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1. Summary of methods. Let A be a non-defective (see §2) square matrix of order n, 
symmetric or not, for which it is desired to determine some of the characteristic values 
vy and associated column-vectors X and row-vectors Y. In terms of matrix products 
these quantities are defined by the relations 


AX = rX, YA = pY. (1) 


A class of numerical procedures is based on iteration methods to obtain one character- 
istic value \ and the associated vectors C, R. Then A is transformed into a matrix A’ 
and a new iteration is used to obtain a characteristic value v’ and characteristic vectors 
X’, Y’ of A’, which can then be converted into corresponding quantities vy, X, Y for A. 
If more values are wanted, one can continue by transforming A’ to A”, and so on. 
Vector iteration schemes for getting one characteristic value of a matrix were described 
in 1929 in [15]; these methods are explained and extended in [1, 13, 8, 9, 5]. 

In the present paper we are not interested in the iteration procedure as such, but 
wish to discuss a class of transformations whereby A’ is obtained from A. The earliest 
of these known to us is “deflation,’”’ suggested by Hotelling [6, 7] for symmetric matrices 
and extended in Aitken’s thorough study [1] to non-symmetric matrices, defective or 
not. In [3] and [4, p. 143] Duncan and Collar introduced a different transformation 
(see §3) for non-defective matrices; this was restated in [10] and [16]. It has the ad- 
vantage that it reduces the order of the matrix, but it destroys the symmetry. In a 
relatively inaccessible paper [13] Semendiaev gave a careful exposition of Aitken’s 
techniques, and extended them tc cover the case of multiple characteristic values in 
full generality. His transformation is very general; in the simplest case it somewhat 
resembles that of Duncan and Collar. Semendiaev expressed his transformation in the 
form of a matrix relation A’ = UAU™. Blanch has devised (unpublished) another 
modification of the Duncan-Collar reduction in the form UAU™. 

In [14] Tucker published a related transformation yielding a matrix A’ of order 
n + 1 which is defective with respect to a double characteristic value zero. Although 
the coefficients are obtained easily by bordering A, the increased order may be a dis- 
advantage. Tucker’s method is not directly a special case of our (5). 

A distantly related matrix transformation is the ‘‘escalator method” of Morris and 
Head [12, 11]. It relates the complete set of characteristic values and vectors of A to 
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the complete set for a submatrix of order n — 1. For this reason the escalator method 
cannot be compared to the transformations considered in this paper, in which at each 
stage one deals with only two characteristic values of A. 

In §2 we present in formulas (5) a four-parameter family of transformations from 
A to A’. This family is general enough to include deflation and the procedures of Duncan 
and Collar, Semendiaeyv, and Blanch as special cases; see §3. In §4 two subclasses of 
the transformation are discussed: order-reducing and symmetry-preserving transforma- 
tions. Two new methods which are both order-reducing and symmetry-preserving appear 
promising for practical work with symmetric matrices A. Even for non-symmetric 
matrices we feel that our family of transformations offers a choice of procedures which 
may occasionally prove useful. 

2. The general reduction method. The reduction formulas will be proved in all cases 
by means of the following lemma, which can be easily verified. 

LemMa. Let the matrix A have the characteristic value v with corresponding column-vector 


X and row-vector Y. Let U be a non-singular matrix. Then the matrix A’ = UAU™~ has 


the characteristic value v with corresponding vectors 
: we : ae 
A’ mm UX. Y moe) (2) 


We assume for simplicity that A is not defective.* Let \ be a (known) characteristic 
value of A, with a corresponding column-vector C = {c, , -:+ , ¢,} and row-vector 
R = (r,,°+:*,7,). Let v be some other (unknown) characteristic value of A, with corre- 
sponding column X = {2,,-:- ,2,} and row Y = (y,, °-: , yn). Assume C, R, X, Y 
to be so normalized that 

RC = YX = 1. (3) 
The characteristic values \ and » are permitted to be equal, provided that X, Y satisfy 
the orthogonality conditions 


C= Rx = 6. (4) 


which are automatic when \ # »v. 
Let y, p, 8, t be complex parameters, and put for abbreviation 


rPel—*¥*,, P=l-p,, f=TP-— vp. 


We now introduce a new matrix A’ = A’(vy, p, 8; t), defined for all values of the param- 
eters except when 8 = 0 and at the same time f + 0: 
at, = Ai; — YCiAnj — PV jdin + VeCT; (Ann — A) — fei; , 
al, = —Bldin — Y¥Cidnn — tev,T + (A — dyei], 
| — fs "fa, — prjAnn — te,r;P + (A — Oopr; (if 8 # 0) ( (5) 
a= 
0 Gf f = B = 0) (Ii<t:<n-lil<j<n-l}), | 
ar, = f(Qan — learn.) + (A — OU — ff). ) 





*In the terminology of [9] a defective matrix A is one for which no transform PAP~ is a diagonal 
matrix. An equivalent definition is that A has one or more non-linear elementary divisors; see [2]. 











No. 


thod 
ach 


rom 
can 
3 of 
na- 
ear 
tric 
ich 











1951] MATRIX TRANSFORMATIONS FOR OBTAINING CHARACTERISTIC VECTORS 327 


For all values of the parameters a matrix U = U(y, p, 8) will be defined below, with 
the property that . 


A'(y, p, 8; t) = U(A — tCR)U™. (6) 


By the lemma A’ has the characteristic values of A — tCR, i.e., those of A except that 
the single characteristic value \ is changed to \ — t. In particular, »v is a characteristic 
value of A’, and the vectors X, Y are transformed into characteristic vectors X’ = 
UX and Y’ = YU of A’ corresponding to v. By use of the formulas for U given below, 


it is easily shown in each case that 


xi = 2 — Wit i<i<n-l), 


(—fx,8"' (if 8 ¥ 0) 
x (7) 
0 af f = 6 = 0), 


Yi = Yi — Yn i sjsa— 0, y, = —By, J 


The vectors C, R are transformed into characteristic vectors C’ = UC, R’ = RU™ of 
A’ corresponding to \ — t. The formulas for the components of C’, R’ are given for 
each case below, as are the definitions of U, U~’. 

An arbitrary choice of the parameters y, p, 8, t may be used to analyze the matrix 
A. Knowing \, C, R, one computes A’ from (5). By iteration or otherwise one next 
determines a new characteristic value of A’. By the lemma » is also a characteristic 
value of A, and the vectors of A corresponding to vy may be calculated from the relations 
(7); these give x, and y, in terms of x{ and y{ except when f = 0. (The case 8 = 0 is 
included in the exceptional case f = 0, as was stated before (5).) When f = 0, one first 
gets x, , y, from formulas (8), which are derived from (7), (4), and (3): 


n-1 
—(y + 7,7)" ) x "2: , 


I 


(8) 


n—-1 


—(o +¢,P)" Do yi: - 


i=l 


Yn 


If still more characteristic values of A are desired, one can use v, X’, Y’ in (5) to trans- 
form A’ into A”, and so on. The procedure is useful to the extent that it is easier to 
find v, X’, Y’ from A’ than it is to find v, X, Y from A. 

It remains only to exhibit U = U (vy, p, 8), so that the reader may verify equation (5). 
For completeness C’ and R’ are also given. There will be three cases. All matrices are 
exhibited in a partitioned form, with a square matrix of (n — 1)-th order at the upper 


left. 
Case 1. T ~ 0, 8B ¥ 0. Here we define 


U = Uy, p, 8) = [errr Fe teeeeeee ees 
—prjB' + —(f + pr,)B™ 
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where 6;; = 0 (¢ ¥ j) and 6,;; = 1 (¢ = j). It can be verified that 
5; — ye! ‘er, —BI ye; 
Pe Ge aa cccduueuenndiue 
— pr;I —pr' 


It is found that 


ci = Tel <‘4<n-— JD), c, = —pB"' — fe,B8'; 


r= fl 'ri(1 < jsn-— 1), ri = —ByI" — Br, . 


Case 2. T = 0. 8 ¥ 0. Here y = c,', and we define 


ae A RE ssa os nner ntaanonssaddsncadeissaneceneys 


6;; — Ci; + pe, (7c, — Dew Bcc, '(r,C, — 1) 
ieee OE esses cite a Pad bin ace Wore en iw Soroe ke aeRO SUES 
—(2p + ¢,P)r —B8 + Bt, — 1)/ 
It is found that 
c=c(1l <i<n — J), c= —(pt+ c,P)B’, 
r} = —pric, (1 < jsn- 1), r, = —Bc;’. 
Case 3. f = 6B = 0. Here y, p are restricted to such values that f = 
(1 — ye,)(1 — pr.) — yp = 0. We define 
6;; — Te, —yc; — Ter, 
U = Uy, p, 0) = rere eee eT eee eT Pee eee . 
r; S 


It is found that 


I 


c=0 (1<i<n-1),, c=], 


r; = 0 (<j<n-— DD), ; = {, 


n 


3. Known special cases of our transformation. For certain values of the parameters 
y, p, 8, t the matrices A’ defined by (5) have previously been used to analyze matrices 
A. We know of the following special cases: 


(a) Duncan and Collar [3] and [4, p. 143]: Case 2, with y = c,’, p = 0, 8 = —1, 
¢ = i. (Rows and columns have been interchanged in our presentation.) The 


matrix A’ is the result of subtracting from each of the other rows c,* times the 
matrix product of C by the last row of A. 
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(b) Hotelling [6] (deflation): Case 1, with y = p = 0, 8 = —1, ¢ = X. Here U is 
the identity matrix. 

(c) Semendiaev [13, p. 212]: Case 3, with y = c,', p = 8 = t = O. This is only a 
special case of Semendiaev’s general reduction. 

(d) Blanch (unpublished procedure used at the National Bureau of Standards, Los 
Angeles): Case 1, with y = 0, p = r,', 8 = 1, t = 0. To compare method (d) 
with (a), to which it is closely related, one should interchange rows and columns. 

4. New special cases of our transformation. One useful class of matrix transforma- 

tions consists of those in §2 for which f = 0; here 6 and ¢ remain unrestricted. For these 
it is seen that af; = 0(j = 1, 2, --- , n — 1), so that A’ is essentially reduced to order 
n — 1. We call these transformations order-reducing; by their use subsequent iterations 
become shorter. The methods (a), (c), (d) in §3 are order-reducing. 

Another special class of matrix transformations consists of those in §2 for which 

y = p, 8 =f, with t unrestricted. When A is symmetric it is reasonable to pick c; = r, , 
and then A’ is also symmetric; hence this class of transformations is called symmetry- 
preserving. In §3 only method (b) is symmetry-preserving. 

New transformations which are both symmetry-preserving and order-reducing are 

those in Case 3 of §2 for which p = y. Except for the unessential freedom allowed ¢, 
there are commonly two of these transformations. When c; = r,; these may be defined by 


c+)" =((%+ 1", (9) 
(¢, — 1)”. (10) 


Symmetric matrices are more convenient to deal with than non-symmetric ones, in that 
by their use the storage requirement is approximately halved and the round-off errors 
are more easily estimated. For dealing with symmetric matrices A, therefore, the 
transformations defined by (9) and (10) look promising. For non-symmetric matrices 
A the method (d) of §3 seems quite satisfactory, but it may occasionally be useful to 
have other subcases of (5) available. 

5. Numerical example. From [8, p. 327] we obtain the symmetric matrix* 


—-2 -—2 0 3 -!l 


zr ? 


= o> (c, ne es 


—2 0 -3 5 0 


3 5 1-3 -1l 








—] 0 1-1 —-1 


By an iteration one can obtain the dominant characteristic value \ = —9.88649 and 
corresponding normalized row-vector 
R = (—.35616, —.52348, —.46374, .61437, .08124). 


Since A is symmetric, the column-vector C has the same components. 





*We have corrected a misprint in [8]. Mr. William Paine of the National Bureau of Standards, 


Los Angeles, assisted with the calculations. 
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We shall apply the reduction defined by (9) to obtain the subdominant characteristic 
value of A. One finds that y = p = (1.08124)~" = .924864. Putting these values of p, y 
into (5) with ¢ = 0, 8 = 0, one finds that 


al; = a;; — .924864(c,a,; + 7;a;,) + 7.60127c,7,, a, =a, = 8 


(lsti<n-1lsjsn-—)), 


a’, = —9.88649. 
Hence 
— 1.69458 —1.06695 1.15597 1.57555 0 | 
~1,06695 2.08298 — .67057 2.07121 0 | 
| 
| | 
A'=| 1.15597 — .67057 —2.50753 —2.16277 0 ; 
1.57555 2.07121 —2.16277 1.00552 0 | 
0 0 0 0 ~ 9.88649 | 


By a separate iteration of the non-trivial fourth-order minor of A’ one can determine 
that vy = —4.75772 is a characteristic value of A’, with corresponding row-vector 


Y’ = (—.54888, —.15891, .67842, .46175, 0). 


By the lemma of §2, —4.75772 is a characteristic value of A, and it remains only to 
determine the corresponding row-vector Y. By (8), 


n-1 
Yn = — Dd, ye; = —.24775. (11) 


t=] 
From (11) and (7) one can then compute the first four components of Y. One finds that 
Y = (—.46727, —.03896, .78468, .32098, —.24775). 


If more characteristic values and vectors of A were desired, one would start by using 
y and the first four components of Y’ to reduce the first four rows and columns of A’ to 
a matrix A” of order four which would be bordered with zeros. 
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1. Introduction. Recently some examples of two-dimensional flow patterns with a 
free stream Mach number one were given. (Refs. 1 and 2) The basis for their determina- 
tion was an investigation of the flow at infinity in the physical plane. It is practical to 
carry out this investigation in the hodograph plane, but in principle the same question 
can also be dealt with in the physical plane. In this form the discussion can be extended 
to the case of axial-symmetry. This investigation is carried out in the present paper. 

2. The flow pattern at Mach number one. Before we proceed to the mathematical 
discussion let us first describe qualitatively the flow pattern whose behavior at infinity 
we are trying to find. Assume we have an axial-symmetric body in a wind tunnel with 
a circular cross-section. If the tunnel is operating under choking conditions, then up- 
stream of the body the flow is subsonic while downstream it is supersonic. These two 
regions are separated by a sonic line which proceeds from the surface of the body to 
the wind tunnel wall. From the supersonic portion of the body two classes of Mach 
waves will arise; those of the first class start from the body and end at the sonic line, 
while the waves of the second also start from the body but end at the tunnel wall. 
There is one Mach wave which separates these two classes, namely, that one which 
reaches the sonic line at its intersection with the wall. We shall label this wave the 
limiting Mach wave. It is clear that the subsonic portion of the flow is not completely 
determined by specifying only the subsonic part of the body. The portion of the super- 
sonic contour along which the Mach waves of the first class originate is also necessary 
since perturbations propagating along these waves will reach the sonic line and thereby 
influence the subsonic region. Therefore, the limiting Mach wave forms the boundary 
of that portion of the supersonic region which is necessary to determine the subsonic 
field. The flow further downstream represents a purely supersonic problem which may 
be computed independently by the method of characteristics. 

If now the diameter of the test section is increased while the body size is kept con- 
stant, then the free stream Mach number upstream of the body will tend to one. There- 
fore, it is plausible that a limiting flow pattern for an infinite diameter of the tunnel 
should exist representing the steady flow with the free stream Mach number one. Here 
also a limiting Mach wave will be present, so that for the determination of the problem 
one need consider only the flow upstream of the limiting Mach wave. 

3. Basic equations. Let x and y be the Cartesian coordinates in a representative 
meridian plane of an axially symmetric flow, the x-axis being in the free stream direction, 
a* the sonic velocity, @ the velocity potential, and x the ratio of the specific heats. 
If one defines a perturbation potential ¢ which gives the deviation from a sonic flow by 


g(x, y) = o(x, y) — a*z, (1) 


then with the approximations usually made in the derivation of the transonic similarity 
law one obtains the non-linear equation 


‘Received July 15, 1949. This is a condensed version of USAF Tech. Report No. 5783. 








G. GUDERLEY AND H. YOSHIHARA [Vol. VIII, No. 4 


wv) 
w 


/ =§ " 
ts + 1) gr¢22 + Guy = ~Y G- (2) 
The behavior at infinity of a flow with the free stream Mach number one will be 


represented by a particular solution of the form 
(zy =y" SO (3) 


where 


' 1/3... —n 


C= («¢ + 1) xy (4) 


and » is a parameter to be determined. The desired solution may be called the basic 
singularity. Inserting this hypothesis into Eq. (2) one finds that f must fulfill the differ- 
ential equation, 


(f? — n°e*)f"’ + (5n — 4)nef’ — (3n — 2)’f = 0 (5) 


Here the primes denote differentiation with respect to ¢. 

The discussion of the integral curves of this non-linear equation is facilitated by a 
group property which it possesses. If a scale transformation is introduced to both f 
and ¢ in Eq. (5), it is found that the differential equation is invariant with a suitable 
choice of the scale factors. For a differential equation with this property it is possible 
to reduce the order of the equation. Indeed if two new variables s and ¢ are introduced 


by 


then one will obtain the first order equation 


dt 2t? + (3n — 4)nt — (3n — 2)*s = 
ie =a 7) 


ds ae (n* — f)i(t— 3s) y 
If one now solves Eq. (7), then one can obtain the corresponding integral of the 
second order equation (5) by means of Eqs. (6). In this connection the value of ¢ is 
found from the integral of Eq. (7) by the following auxiliary equation derived from 
Eqs. (6), namely, 


d log ¢ _ d - (8) 


et 8 6t=—8 
Actually the s, ¢ representation will be used solely for the discussion of the integral 
curves while the final computation is carried out by means of a numerical integration 
of the differential equation for f. 

4. Structure of the integral curves in the s,/-system. A survey of the integrals of 
Eq. (7) is given in Fig. 1. For this purpose the entire s,t-plane has been projected onto 
a hemisphere as shown in Fig. 2 and the hemisphere is then distorted to form the plane 
circular region of Fig. 1. Infinity in the s,f-plane is then represented by the outer circle. 
The arrows on the integral curves show the direction in which the absolute value of £ 


increases. 


As is seen from Eqs. (3) and (6) the sonic line (g, = Q) is represented by the line 
t = 0, that is, line B’AB in Fig. 1. The subsonic region lies below this line while the 


supersonic region lies above it. 
Along the line ¢ = n° one finds an extremum of ¢. This may be seen by inserting 











(2) 


be 


le 
le 
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the differential equation (7) into Eq. (8) and then setting ¢ = n’. Let us examine the 
significance of an extremal value of ¢. If one travels along an integral curve of the 
s,t-system, then the corresponding lines ¢ = constant sweep a certain region of the 
physical plane. The occurence of an extremum of ¢ along the integral curve will indicate 
that one sweeps out a certain region of the flow pattern up to the extremal value of 
¢ and then resweeps the same portion of the physical plane. One thus obtains a multi- 


valued physical plane, the line for which ¢ is an extremum being a limiting line. Conse- 








Fic. 1. The integral curves of the s,f-system projected upon a hemisphere for n = 4/7 


quently all integral curves which cross the line ¢ = n° (dashed line in Fig. 1) are of no 
physical interest. As will be seen later, the integral curve which goes through point C 


will be an exception. 
Let us consider the singularities of the differential equation (7). The coordinates of 


the singular points in the finite portion of the s,t-plane are 


n*(5n — 4)(3n — 2)’, s = 2/9, 


s = 0, Ss 
t = 0, t=n’, t= 2/3. 


They are shown in Fig. 1 at points A, C, and D, respectively. 

Point A is a nodal point for all values of x while points C and D may be either a 
nodal point or a saddle point depending upon the value of n. For our problem n will be 
confined to the range between zero and one because of the requirement of zero pertur- 
bation velocities at infinity. For 0 < n < 2/3 point C will be a saddle point and point 
D will be a nodal point while for 2/3 < n < 1 the two points interchange their structures. 
For n = 2/3 the two points coincide to form a singular point which has both the char- 
acter of a nodal point and a saddle point. 

The singular points at infinity are investigated best by means of the introduction 
of homogeneous coordinates, 2) , 2, , and x, of Fig. 2. One finds two singular points at 
infinity. Point B (which is identical with B’) is a combination of a saddle point and a 
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nodal point; it is reached as ¢ tends to infinity at a constant value of s. Point FE is a saddle 
point, and it is reached by traveling outward along a line s/t = 1. 

The significance of the different singular points is seen if one integrates the differential 
equation approximately for the neighborhood of the singularities and then determines 
the behavior of ¢ and the velocity components of the perturbation potential 

It is found that ¢ tends to infinity as one travels towards point A along any of the 
integral curves. Therefore point A corresponds to the x-axis. Among the integral curves 
which enter point A only one can be expressed in terms of a power series. All other 
integral curves will contain a logarithmic term. It is found that only for this one in- 
tegral curve will the velocity component in the y-direction tend to zero as one approaches 
point A. 


























Fic. 2. Orientation of the homogeneous coordinate system. 


For all integral curves which enter point B ¢ tends to zero; that is, point B corre- 
sponds to the y-axis. The velocity components at this point will in general assume finite 
values with no special restrictions imposed upon them. Since a subsonic velocity cannot 
jump suddenly to a supersonic velocity or vice-versa, the continuation of a curve which 
enters point B with a subsonic velocity, for example, must leave point B’ again with 
subsonic velocities. 

At point C ¢ assumes a finite value. If one rewrites the expressions for the coordinates 
s and ¢ at this point in terms of the components of the perturbation velocity, then one 
finds that point C corresponds to a line ¢ = constant which coincides with a Mach 
wave. Furthermore, it can be shown that all Mach waves which lie upstream of that 
represented by point C will reach the sonic line in a finite distance from the axis. There- 
fore, point C will represent the limiting Mach wave described in Section 2. 

The other singular points are without interest in the present case. 
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5. The determination of the paramenter n. It is now possible to describe quali- 
tatively the properties of the integral curve which will represent the desired solution. 
The integral curve will start at point A which represents the horizontal axis upstream 
of the body. Only that solution which is regular here can be used since along the negative 
z-axis the y-component of the velocity is obviously zero. One will then proceed in the 
direction of subsonic velocities to point B which corresponds to the y-axis. The in- 
tegral curve now reappears from point B’ with a subsonic velocity and then crosses 
the sonic line B’A and enters the supersonic region. It should end at point C which 
represents a line ¢ = constant corresponding to the limiting Mach wave. Naturally for 
this purpose it is necessary that the value of n be chosen correctly. 

It is interesting to see what happens to the solution in the case that the value of n 
has been chosen incorrectly. The course of the integral curve from point A to point B 
will be qualitatively the same for values of n in the neighborhood of the desired one. 





Fic. 3. The functions f, f’ and f” 


The main symptoms of a wrong choice of » will appear in the continuation from point 
B’. Depending upon the choice of n one will have curves either above or below the de- 
sired curve B’C. The integral curves which lie above B’C will intersect the dotted line 
of Fig. 1; they will therefore represent flows with a limiting line. Those integral curves 
which lie below B’C will cross the sonic line B’A and enter point A with supersonic 
velocities. One will expect that in this case one does not enter point A along that in- 
tegral curve representing zero downwash velocities along the x-axis. 

Are such flow patterns physically possible? Let us consider the first case in which a 
limiting line occurs. In this case, the limiting line is formed by an envelope of expansion 
waves. Such a configuration of Mach waves in which one expansion wave starts from a 
neighboring one in the interior of the flow may be considered as arising from a series of 
thin knife edges placed along the limiting line (see Ref. 3). A flow pattern of this kind 
arises for a value of n for which the perturbation velocities vanish rather rapidly as one 
goes to infinity. 

The integral curve which lies below the desired one, on the other hand, represents a 
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flow in which the perturbation velocities decrease very slowly towards infinity. Here 
the flow will correspond to that over a semi-infinite body whose contour is shaped in 
such a way that the expansion waves originating from the body are never sufficient to 
produce supersonic velocities up to a very great distance from the body. In this case 
the limiting Mach wave does not lie in the finite portion of the physical plane, so that 
all Mach waves from the body will reach the sonic line. 

According to this description the flow patterns of the first type have a high density 
of rarefaction waves which can be produced only by a limiting line, while the flow pat- 
terns of the second have a sparse distribution which requires a carefully shaped body 
which extends to infinity. The desired flow represents the transition between the above 
two types of flows of opposite character. 

As mentioned previously it was found more practical to carry out the numerical 
computation in terms of the differential equation for f. The most satisfactory procedure 
for the determination of n was to begin the integration from both points A and C and 
then match the two resulting integrations at a point intermediate to the two points 
using the s,f-representation. The solutions near the singular points are found by a suit- 
able series representation. Within the accuracy of the numerical computation the value 
of n was found to be 4/7. The resulting function f and its derivatives f’ and {” are shown 
in Fig. 3. The streamline pattern for the flow represented by the basic singularity is 


MY AQ 

es oF | _| 7 : as | 

SVT 
ay 


saga \ LLY 
AAT 
\ VI Aspe CONSTANT 


_ V /} Py 
NA / 


} 
-4 -3 -2 — fe) ! < 3 4 


shown in Fig. 4. 


































































































Fic. 4. Lines of constant @ and the streamlines for the flow represented by the basic singularity. 


6. The final results. The basic singularity given by Eq. 3 with n = 4/7 may now 
be used to study the perturbation velocities at infinity caused by an axial-symmetric 
body with a free stream Mach number one. (Here the role of the basic singularity is 


analogous to that of a doublet singularity in the case of a body in an incompressible 


flow.) 
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The results show that as one proceeds farther and farther from the axis of sym- 
metry along lines of constant ¢, where ¢ = (« + 1)~’*xy"*”’ the perturbation velocities 
respectively. (In the two-dimensional case these 


g- and y, decrease as y“’’ and y*”’ 
velocities decreased as y*’’ and y"*’’ as one traveled outward along the curve zy~*” = 


—6/7 


constant. Cf. Ref. 1.) 
In Fig. 5 is shown the distribution of the pressure coefficient C, and the streamline 
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deflection @ along a horizontal line at different distances from the axis. Here C, = 
2(p — p*)/p*a**, where p and p are the local pressure and density and the asterisk 
denotes the sonic free stream conditions. 

It is of further interest to determine a body shape which produces a flow at infinity 


given by the basic singularity. An approximate construction of such a body has been 
carried out in Ref. 4, and the resulting contour is shown in Fig. 5. 
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THE RHEOLOGICAL ASPECT OF HYDRODYNAMICS* 
BY 
M. REINER 
Hebrew Technical College, Haifa 
1. When the rheologist considers liquids, he is interested in them as materials having 

certain material properties. Now, when he turns in his investigations to hydrodynamics 
for help, what she offers is the mathematics of two sorts of liquids:—One is called the 
ideal liquid and has no material property besides elasticity of volume determined by the 
bulk modulus «, connecting a hydrostatic pressure p with a volume expansion e, by 

p= —Ke,; (1) 
the other, the viscous liquid, has, in addition to volume-elasticity, a property called 
viscosity (u),’ defined quantitatively since Newton as the ratio of shearing, or tangential, 
stress (p,) to the velocity gradient (@), or 


k= pr G. (2) 


The latter is better replaced by the tangential component of the rate of strain (e,) which 


is half the velocity gradient” 


e, = G/2; (3) 
so that 
e, = p,/Qu. (4) 


The rheologist is at once confronted with the question: Is this enough? Can he for the 
description of the rheological behaviour of liquids manage with these two “coefficients”: 
« Which refers to afl materials; and uw which refers specifically to liquids? This question will 
presently be examined. 

2. Rheologists very often have to deal with highly viscous substances, of the order 
10° to 10°° poises, sueh as tar and pitch, for which the standard “shear’’ methods of 
viscometry are very inconvenient. Since Trouton (1906) such substances are formed 
into cylindrical rods which are subject to pull and the rate of elongation is observed. 
Let p.. be the tractional force per unit area of cross section and e,. the rate of elongation 
per unit length, with the z-coordinate in the direction of the axis of the cylinder, then 
Trouton defined a coefficient of viscous traction by the equation 


Ar = Pis/ere (5) 


The coefficient \,7 will evidently depend in some manner upon uz. As a matter of 
fact the Stokes-Navier differential equations, which are claimed by classical hydro- 
dynamics to be all-embracing, should be sufficient to deal with the case. 


*Received May 6, 1949. Lecture given on October 10th, 1946, at the Cavendish Laboratory in a 
colloquium arranged jointly by the Department of Colloid Science, Cambridge University, and the 
Cavendish Laboratory. 

‘For easy reference I shall use Lamb’s notation. Rheologists usually denote the coefficient of viscosity 


by ». 
2e’, is sometimes identified with G, but (3) is a better definition. 
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In order to find the relation between \7 and » Trouton reasoned as follows:* 

When the pull p,, is applied, there will be an /nitial state, the first stage in the ex- 
periment, when the rod is not only elongated but at the same time also expands, the 
measure of the volume expansion being at every moment e, and its rate e, = 6. This 
volume expansion produces an elastic reaction of the nature of a hydrostatic tension 
in accordance with (1). The pull p., has an isotropic component p,, 3.* When e, has so 
much increased that xe, = p,,/3, the elastic reaction balances the isotropic component 
of the pull. With this, the volume expansion ceases to increase and its rate of increase 
e, accordingly vanishes. A second stage then sets in, in which e, = Ois permanently equal 
to zero. Let us forget about the first stage and start observations only after the second 
stage has been reached. We then can make use of the equations (4,L) Art. 328 on page 
577 of the sixth edition of the classical work of Lamb (1932). In order to get rid of the 
action of gravity, we immerse the rod, as was done by Trouton, in a liquid of equal 
density. The movement being slow, we may also neglect inertia forces. Then the Stokes- 
Navier equations become 

—dOp/dx + uVou = 0 (6) 
and two others where z is replaced by y and z respectively and u by v and w. 
A solution of (6) is given by 
u= Qn, v= qy, w= 1r2Z (7) 


in which case p, the hydrostatic pressure, is constant throughout the rod. 
In order to determine g and r, we go back to Lamb’s equations (5, L), and (6, L), 
Art. 326 on page 574, from which we find the stresses 





Dez = —p + 2pdu/dz; py, = w(dw/dy + dv/dz) (8) 
and similarly p,, , Pez» Dez » Dev + 
Now from our (7), 
du/dx = dv/dy = g; Ow/dz = 7; dw/ dy + dv/dz = D0 (9) 
and therefore 
2 
Di = Dey = “PP — ST AO Ds 
ow 
(10) 
Dee = —Pt 7 y — @). 


We know p,, and as no forces act on the sides of the bar, p,, and p,, must vanish. This 


gives us two equations 


0) 
2u 
p+ 3 (r — g) = 0, 
(11) 
_— tu lr — = _— 
p = v qg) = ae 
vw 
In what follows Trouton’s arguments are not reproduced literally, but in their essence. 
‘The isotropic component of a stress is the mean normal traction pm = (pPrz + Pyy + Pzz)/3. In our 


Case Prr = Py = 0. 
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which together with the third 


c, = 6 = du/dx + dv/dy + dw/dz = 2g +r =0 (12) 
yield 
p = —p,,/3 (13) 
and 
r= p,,./3u, 
(14) 
q = —D./6p, 
or [compare (9)] 
C.. = Ow/dz = r = p,,/3u (15) 


from which [compare (5)] 
Ay = 3p. (16) 


3. However, this reasoning of Trouton’s begs the problem. How is one to know when 
the first stage has come to an end? In materials of such high viscosity this may, for all 
we know, take quite a long time. One cannot rest satisfied with such a state. The rheol- 
ogist often has a practical interest in the first stage. He has the feeling that Stokes 
conjured away another coefficient of viscosity. Actually Lamb’s equations (5, L) and 
(6, L) which we used above, are gained by specialization of more general equations, 
namely. (2, L) and (3, L), page 574 which are® 


Pir = —DP +O + Qy du/dz, 


(17) 
Dye = w(Ow/dy + dv/dz). 
By defining (compare Art. 325, p. 573) 
—P = (Dez + Pw + Des)/3 (18) 
and introducing this definition, Lamb finds 
3\ + Qu = 0 (19) 
from which 
A= —2y/3. (20) 


This must cause confusion. As p,, , P,y » Pzs May, on principle, depend “also on the 
rate of expanson at the point (z, y, z)°”’, p, accordingly, is in general not the hydrostatic 
pressure, i.e. not that pressure which is independent of the rate of strain. If we, as most 
workers will be inclined to assume without question, use p to denote the hydrostatic 
pressure, we shall not find (19), but will in general have 


3\ + Qu = 3y,, (21) 


’This \ should not be confused with Trouton’s coefficient of viscous traction Ar . 
‘Compare footnote on page 578 1.c. 
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where uw, is another coefficient which may be named ‘‘volume viscosity’’. It therefore 
appears that Stokes had assumed the volume viscosity to vanish. This is borne out 
by the following quotation from Stokes, to which Tisza (1942) has drawn attention: 

. of course we may at once put uw, = 0, if we assume that in the case of a uniform 
motion of dilatation the pressure at any instant depends only on the actual density 
and temperature at that instant and not on the rate at which the former changes with 
the time. In most cases to which it would be interesting to apply the theory of the 
friction of fluids, the density of the fluid is either constant or may without sensible 
error be regarded as constant, or else changes slowly with time. In the first two cases, 
the results would be the same and in the third nearly the same, whether yu, were equal 
to zero or not. Consequently, if theory and experiments should in such cases agree, 
the experiments must not be regarded as confirming that part of the theory which relates 
to supposing uw, to be equal to zero’. Tisza has objected to this assumption, which, 
according to him, is valid in monoatomic gases only while in polyatomic gases and 
liquids u, > yw. His arguments are based on considerations of supersonic absorption. 
We shall presently show that at the other extreme in very viscous liquids Stokes as- 
sumption also leads to unacceptable results. 

4. We now start from 


Diz = —D + (u, — 2u/3)0 + Qpdu/dr (21) 


Only p,,. , the ‘viscous traction” acts in our case and both p,, and p,, are absent 
Therefore 


—p + (u, — 2u/3)6 + Qu du/dx = 0, 
—p + (u, — 2u/3)6 + Qu dv/dy = 0, (22) 


—p + (u, — 2u/3)60 + Qu dw/dz = p, 
from which, by adding up, 
6 = (p../3 + p)/bM - (23) 
Introducing this expression into the third of (22) we get 
e.. = 0w/dz = [p ee — ig, = oe 3)0|/2Qyu (24) 
which on re-arrangement becomes 


3p :2h,|/Gube - (25) 


e.. = [(p.. + 3p)yu + 

Introducing the expression for the hydrostatic pressure p from (1) we now get for (5) 
Ar = e,,/P.2 = [(1 — 3xe,/p..)m + 3u.]/9un, (26) 

and therefore 

Ar = [guu.(1 — 3xe,/p 2 gt fu + 3u,(1 — 3xe, '»..)'). (27) 
Equation (27) is the viscous analogy to the classical relation between Young’s modulus 
(E), the modulus of rigidity (G@) and the bulk modulus (x), viz. E = 9G@x/(G@ + 3k). 
It should, however, be noted that the analogy is not exact. While FE = 3G for x = 
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and only in this case, the situation is very different with regard to \,7 as will presently 
be discussed. 

5. In interpreting (27), one should keep in mind that 1 — 3xe,/p., 2 0, where the 
sign of inequality is valid in the first, the sign of equality in the second stage of Trouton’s 
experiment. For the volume viscosity we have 0 S yp, S@. It is clear that the volume 
viscosity cannot be negative as ‘‘otherwise the more alternate expansion and compression, 
alike in all directions, of a fluid, instead of demanding the exertion of work upon it, 
would cause to give work out” (Stokes). As can be seen, in the first stage, if u, vanishes, 
\, also vanishes. This would mean that in the first stage of Trouton’s experiment the 
viscous resistance of a very viscous liquid against extension vanishes, no matter how 
high the ordinary viscosity w of the liquid, a result at variance with our ideas of viscous 
flow. On the other hand \;/u= 3 in both stages if 4, = © and also in the second stage 
whatever the magnitude of yu, . Because u, cannot be negative, \, cannot exceed the 
value 3. Therefore 0 S Ar S 3. 

Examination of (27) now shows that Stokes was mistaken in equalizing the influence 


of either e, = 0 or wu, = 0 on experimental results. We see that the same result follows 
from either e¢, = 0 or un, = © and not uw, = 0 and generally u, > uw. It may be mentioned 


that Tisza deduces from supersonic absorption in certain liquids for u,/u a value of 
2000. This, for all rheological purposes, is infinite. 

6. Having discovered in this way arother rheological coefficient besides u, namely 
u, , we have become suspicious. There may be many other coefficients, concealed from 
us by the biased procedure followed in classical hydrodynamics. The rheologist therefore 
makes himself independent and starts ab ovo. 

We first have to take account of the fact that even the simplest liquid is complex 
insofar as it has something in common with the solid, namely the property of not flowing 
under the action of a hydrostatic pressure (or tension). A hydrostatic stress will change 
the volume of the liquid in a reversible elastic manner. That stress is equal to xe, . In 
order to investigate the rheological properties of the liquid qua liquid, we must deduct 
from the applied stress p*,’ the hydrostatic stress xe,6,, where the stress tensor p,, is 


defined by 


Dez Dey Des | 
|| 
Pes = || Pur Pov ~~ Pv || (28) 
| 
| Dez Prey Pez | 
and the unit tensor 6,, by 
li 0 0 
| 
5.=|10 1 O (29) 
| 
| | 
0 0 || 
Only what remains, namely 
Pro = pe, — Ke, 5, (30) 


7The asterisk indicates the stress as applied, i.e. including the hydrostatic part. 
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will be the stress due to flow and depending upon the tensor of rate of strain e,, . Thus 
p,. is the viscous resistance in the widest sense. Instead of the very special equation (4) 
we now write 


Ce = fpr) (31) 


and our task is to develop the function f in as general a way as possible. Reiner (1945) 
has shown how this can be done and the result is 


Cre = Fo(P)6,, + F.(P)p,. + Fo(P)p,aPas (32) 


where the ¥ are functions of the three invariants P, , P, , P; of the stress tensor p,, and 
PrePas is a tensor of the same rank as p,, , the components of which can be written 
down if we remember the summary convention of tensor analysis and note that the 
suffix a appears twice.* 

All rheological coefficients are contained in the furctions $) , $, and 5%, . In order to 
extract them from there, we perform rheological experiments and adapt (32) to the 
experimental conditions. Let us first drop the last term on the right side of (32) and 
deal with the equation 


Cre = Fo(p)5,, + F(P)p,. (33) 


where $, and , are furthermore to be functions of the first and second invariant only. 
It will facilitate our considerations if we resolve e,, , as can be done with every tensor, 
in an isotropic (e,) and a deviatoric (e,’) part as follows 


e, = SP, , P2) (34) 


ns = Fi(P2)p%, , (35) 


Crs 


where P is the second invariant of the deviator, its first invariant being equal to zero. 

We require that the viscous resistance vanishes with vanishing rate of strain and 
vice versa. Equ. (35) conforms to this requirement, but in order to make (34) conform 
to it we re-write it as follows 


e, = P,5i,(P, , Pi) + Pi5i.(P; , P), (36) 


vanishing stress carrying with it vanishing invariants P, and P} . Analysis of (35) and 
(36) will yield us then at least three rheological coefficients. We know up to now two 
(u and u,) and we therefore look forward with interest of meeting a third one. 

7. (i) The first of our experiments, carried out in imagination, will consist in the 
application of an isotropic stress p.* Deducting xe, , we get 


p= p* — xe, (37) 
8 : a=z 
PraPas = =. PraDae = DesDeo 1 DevDus F PrsPeo « 


The first component of this tensor is accordingly 


DezDez + DeyPy2 + DesPez = Dez + Dry + Dis 


and similarly the others. 
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while p/, = 0. Equ. (35) gives e/, = 0 and the rate of strain will have an isotropic com- 


ponent only. The first invariant of the stress pé,, is 


P, = 3p 


(38) 


while the second invariant of the deviator of stress, there being no deviator, is equal to 


zero. This makes (36) 
ce, = 3p%i,(p, 0). 

354, (p, 0) is our first rheological coefficient. We rewrite it in the form 
[355.(p, OT’ = wu, 


and have 
p _ HC» 


and (compare (30))° 
p* = xe, + we 
(ii) For the second experiment we apply a tangential stress p, 
| 0 Dp: 0 | 
| | 


p*, = | Dp: 0 0 
| 


| 





‘0 0 0 
from which (compare (30)) 
ot. Dp: 0 
Paw = De — KC, 0 
0 0 — Ke, 
The first invariant of p,, is 
P, = —3 xe, 
the second of its deviator (compare (43)) 
Pi = —p 
Equ. (35) gives us 
0 Dr 0 


e,= F(p)||\p O 0 


0 0 0 
¥{(p:) is our second rheological coefficient. We re-write it in the form 
[Fi(pi) |" = Qn 


*Such a rheological equation is mentioned by Reiner (1948). 


(39) 


(40) 


(41) 


(42) 


(43) 


(44) 


(45) 


(46) 


(47) 


(48) 
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and have 


DP: = Que, (49) 
in accordance with (4), or generally 
pl, = Que,’ . (50) 
Equation (36) on the other hand gives 
e, = —3xe,5i,(—3xe, , —p) — p?Fo(—Bxe, , —p?). (51) 


This is something new. While in accordance with the assumptions of classical hydro- 
dynamics a simple shearing stress produces a continuously increasing simple shearing 
strain—and nothing else (compare (4)), we now found that generally it will, in addition, 
produce a volume expansion, the rate of which is given by (51). This volume expansion 


raises an elastic resistance which ultimately stops it. When this is the case, we have 


e; = O and 
e, = —(p;/3x)Foo/ For (52) 
This gives us a third coefficient” 
6 = —35/,/F/. (53) 
so that 
= p;/Ké (54) 
or generally” 
= —P; KO (55) 


We have found a third coefficient by mathematical analysis and after we found it, 
we may first be puzzled about what it should signify—but not for long. The phenomenon 
is well known. An increase of volume caused by a shearing stress is what Reynolds first 
described in wet sand as dilatancy. We may, therefore, call 6 the coefficient of dilatancy. 
The phenomenon has more recently been noted in very concentrated suspensoids."* 
Classical hydrodynamics could not account for it. Our analysis has shown that it may 
generally be present in any liquid. 

8. If 5, is a function of the third invariant of stress as well, (36) will become 

e, = P,5i(P, , Pi , Ph) + Pi5i(P, , PsP) + Pisis(P, , PiP?) (56) 

This implies a fourth rheological coefficient 5;, . If we furthermore take up 5, as 
well, we see that a viscous liquid will in general be characterised by 5 coefficients (in 
addition to its elastic bulk modulus) of which, so far, we have defined three only. This 
does, of course, not mean that every viscous liquid will possess rheological properties 
requiring for their description finite values of all five coefficients. Sometimes one and 


WThis 6 is not to be confused with the 6 of 6,, . 

“One should keep in mind that the second invariant of a deviator is always negative. 

12This is the place for a warning against an indiscriminate use of the term dilatancy. It is sometimes 
“structural” viscosity. This is mis- 


identified with an increase of viscosity with shear—the opposite of 


leading. The so-called snow-plough effect in suspensions is not dilatancy. 
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sometimes another property will be more prominent, depending upon the conditions 
under which the liquid is observed. We have had a striking example with the coefficient 
x, - From the times of Stokes to recent times viscous liquids were observed under condi- 
tions which concealed its presence. Only when the rate of strain was very great, as in ultra- 
sonic waves, or the material observed very viscous, as in the case of tar, did it become 
necessary to take this neglected coefficient into account. It may well be that particular 
conditions of observations will require the consideration of other coefficients, to be 
derived by the method shown in this paper. This may be the case with the liquids which 
Weissenberg (1946) has called “general’’ liquids. As was shown by Reiner (1945), the 
presence of F, implies that a simple shearing stress causes not only a continuous shearing 
of the liquid but also a ‘‘continuous lengthening (or shortening) in the direction normal 
to the plane of shear’. Reiner adds that such a phenomenon “has never been observed’’. 
Actually the phenomenon had at that time been observed by Weissenberg and collabo- 
rators but published only later. 

Conclusions. Tensor analysis of the most general relation between stress and rate of 
deformation leads to five coefficients of viscosity of which three are the ordinary (shear) 
viscosity u, the coefficient of volume-viscosity 4, and the coefficient of dilatancy 4, 
while one coefficient may be connected with the ‘Weissenberg-effect’’. The Stokes- 
Navier equations of hydrodynamics in which y, is neglected, can be used consistently 


only when the volume of the flowing liquid is constant. 
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—NOTES— 


QUASI-STATIONARY AIRFOIL THEORY IN SUBSONIC COMPRESSIBLE FLOW* 
By JOHN W. MILES (University of California, Los Angeles) 


Summary. A solution of the integral equation for an oscillating, two-dimensional, 
thin airfoil in a compressible flow (subsonic and inviscid) is obtained by retaining only 
first order terms in frequency. The results are applied to the calculation of the damping 
derivative of a tail in rotary motion about a forward center, and it is shown that the 
damping is considerably less than that calculated on the basis of stationary airfoil 
theory. A brief investigation of induction effects shows this reduction to be considerably 
less for a wing of finite aspect ratio. 

1. Introduction. The problem to be considered in this paper is the determination of 
the pressure distribution on a two-dimensional thin airfoil in a compressible flow for 
the case where the prescribed downwash distribution exhibits a low frequency, harmonic 
time dependence. (By “low frequency”’ we imply that only terms of first order in the 
reduced frequency parameter, k, will be retained.) The situation of primary interest is 
subsonic, compressible flow, since the incompressible’ and supersonic’ cases have been 
considered previously. 

In studying the disturbed motion of an airplane about its equilibrium position, the 
aerodynamic forces generally have been calculated on the basis of steady flow. The 
argument advanced is essentially that the motions are “slow” (i.e., low frequency), so 
that aerodynamic phase shifts may be neglected. The implicit assumption is that only 
first order terms in frequency need by retained and that these terms occur only in the 
(airfoil) displacement. Unfortunately, first order terms also occur in the hydrodynamical 
equations and will be important, insofar as zero order terms occur in the prescribed 
downwash. 

Also, it must be expected that compressibility effects will be rather more important 
in non-stationary flows, since the effect of a finite velocity of propagation on the dis- 
turbances introduced by the airfoil is to emphasize phase shifts, as represented by both 
first and second time derivatives in the potential equation. Thus, the linearized (subsonic) 
compressibility correction of the quasi-stationary results for incompressible flow over a 
thin airfoil is not effected by simply introducing the Prandtl-Glauert factor (1 — M’)~”. 
The variation with Mach number is particularly important at supersonic speeds, where 
the first order (in frequency) terms neglected by steady flow theory may cause a one 
degree of freedom system to pass from an unstable to a stable state (or conversely) at 
the transition of some critical value of M. 

The analysis of quasi-stationary, compressible flow is quite naturally a limiting case 
of the oscillating airfoil theory first considered by Possio.* Fortunately, a closed form 


*Received July 18, 1949. The results in this paper were originally announced in the Readers’ Forum, 
J. Aero. Sci. 16, 511 (1949), and the present material is based on work done at the Aerophysics Lab- 
oratory, North American Aviation, Inc. 

1J. W. Miles, Quasi-stationary thin airfoil theory, J. Aero. Sci. 16, 440 (1949). 

2J. W. Miles, On harmonic motion at supersonic speeds, J. Aero. Sci. 16, 378-379 (1949). 

3C. Possio, L’ Azione aerodinamica sul profilo oscillante in un fluido compressible a velocita iposonora, 
L’Aerotechnica, 18, 441-458 (1938). 
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solution may be given for this limiting case, whereas this has not yet been possible for 
the more general situation. 

2. The Possio integral equation. A thin airfoil is located near z = 0,/2|] < 1 ina 
subsonic, compressible flow of velocity U directed along the positive x axis, as shown 
in Fig. 1. If w(x) exp (wt) is the prescribed downwash at a point x on this airfoil, y(x) 








a i~ hae 
~ Rincon — x 
7 —  % 
@= 0 





3 


Fig. 1. Thin Airfoil in two-dimensional flow. 


exp (7wt) represents the local vorticity associated with this downwash distribution, / 
is the Mach number, and k is the reduced frequency, the required integral equation 


between w and y may be written’* 


aa 
wr = (2r) "(1 — M°*) " | dé (x = ¢)~ exp | —sk(z = &) 
Jn 
(2.la) 
-R[M, k— M2)"e — Bh) 
y(1) = 0 2.1b 
— (wh [ 2 2) 
R(M, y) = : imMy | du\|u\~' exp (tu)H{’(M | u |) 2.3) 


where H;”’ is a Hankel function in Watson’s notation.” Eq. (2.1b) is the Kutta condition 
(continuous pressure at the trailing edge) and is required in order to render the solution 
unique. The semi-chord, }, is the unit of length, as indicated in Eq. (2.2) 

A discussion of Possio’s integral equation and the various methods which have been 
proposed for its solution, together with a discussion of the analogous supersonic problem, 

‘H. G. Kiissner, General airfoil theory, NACA T.M. 979 (June, 1941). 

5G. N. Watson, Bessel-functions, Macmillan Co., New York (1945). 
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has recently been given by Biot and others.° The procedure to be used below is essentially 
> » . 7 . Ps 

that previously proposed for the case of small Mach number.’ In carrying out this 

solution we find it expedient to introduce the modified downwash and vorticity functions 


w(x) = (1 — M’)~'? exp (—iM’*kz)w(2) (2.4) 
x(a) = exp (—iM’*kzx)y(2) (2.5) 
k = (1 — M’)"k , (2.6) 


whence Eq. (2.1) transforms to 


w(x) = (2r)™' [ dé (x — £)~* exp [—ak(x — &)JR[M, k(x — 8] (2.7a) 


y(1) = 0 (2.7b) 


3. Approximate solution to Possio equation. We now set ourselves the problem of 
solving Eq. (2.7) after neglecting terms of order k* or higher. The first step in this solution 
is naturally the expansion of R(., y) in powers of y, with only terms of order y (by order 
y", we also imply terms of order y" In y) to be retained. 

As it stands, R(.M, y) is improper for positive values of y due to the singularity of 
order u~* in the integrand. Actually this difficulty is synthetic and may be removed 
by integration by parts to obtain a proper integral.” Thus, integrating Eq. (2.3) twice 
by parts, evaluating the infinite integral with the aid of the Weber-Shafheitlin integrals 
given by Watson,” multiplying the result through by exp (—7y), and expanding in 
powers of y, We may approximate Eq. (2.7) in the form 


w(x) = (2x) [ dé [(x — &)7' + kh(x — SD] yr) + O}9%) (3.1a) 
y(1) = 0 (3.1b) 
h(x) = 34 — aly + In(k {| 2x |) + f(D) (3.2) 


i(M) = [1 —  — M*)'?] nM 4+ (1 -— M*)'? In (1 + (1 — M’)'?] — In2_ (3.3) 


where y (without functional dependence, to distinguish from y(x)) is Euler’s constant. 
Eq. (3.1) can be solved by iteration, for if y(x) is written yo(2) + Ky,(x) and powers 
of k equated, the results may be identified as the well known integral equation for thin 


6S. N. Karp, 8. S. Shu, H. Weil, M. A. Biot, Aerodynamics of the oscillating airfoil in compressible 





flow, F-TR-1167-ND, HQ, AMC, Wright Field, Dayton, Ohio (1947). 


7J. W. Miles, A note on a solution to Possio’s integral equation for an oscillating airfoil in subsonic 


flow, Q. Appl. Ma. 7, 213-217 (1949). 


8’We do not attempt complete rigor here, but we remark that this difficulty is common to many 
problems of mathematical physics which are formulated as Fredholm integral equations. It usually may 
be removed by a more thorough derivation, including, in particular, the linearization of the boundary 
conditions and the associated passage to the limit after appropriate integrations by parts (rather than 
before, as is implicit in the present case). 

loc. cit. 5, p. 405, Eqs. (4) and (5). 
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airfoils in steady, incompressible flow. While this would be the most direct approach, 
it may be observed that our problem may be reduced to an equivalent incompressible 
problem (neglecting terms of order k’) if k and Euler’s constant (y) are replaced by k 
and y + f(M), respectively. Expanding Schwarz’s solution for the incompressible 
problem,’’ retaining only terms of order k, and substituting k, w, and y, from Eqs. 
(2.4)-(2.6), in the resulting solution to Eq. (3.1), we obtain the desired approximate 
solution to Eq. (2.1) in the form 


ar 


y(— cos #=) = r ‘(1 — MM’)? | dg g(8, g)w(— cos ¢) (3.4) 


g(9, ¢) = 2[cot (6/2) + (cos 6 — cosy)" sin 6] 


f — ‘OS 
+ ik} sin g In }} cos (8 + 2 | 


1 — cos (0 — ¢) 





(3.5) 
+ 2ly + 47m + In (3k) + f(M) — M’] 


° (1 — COS ¢) cot ao a O(k’) 


4. Lift and moment coefficients. The pressure coefficient on the wing is given by 
c,(x) = (4pU*)“"[p(z, O—) — p(x, O+)] = 2U~'y(x) (4.1) 
The lift and (quarter chord) moment coefficients, referred to the chord and the square 


of the chord, respectively, are given by 


c, = U" | d@y(— cos 6) sin 6 (4.2) 
r\ =—1 or . l 
c, = (2U) | d@y(— cos 6) sin 6 cos 6 — 4 Cc; (4.3) 
The downwash due to a displacement z(x), positive down, is given by 
a “oe , 
w(r) = Di 2(x) = U[z’(x) + tke(x)] (4.4) 


Combining Eqs. (3.4), (3.5), and (4.4), substituting in Eqs. (4.2) and (4.3), separating 
the real and imaginary parts, and introducing k from Eq. (2.6), we obtain 


Ri(c;) 


(1 — uM)" | de (1 — cos y)z’(— cos ¢) + O(k) (4.5a) 


Im(c,) 2k(1 — my | dp (1 — cos g)z (— cos ¢) 
0 


+ 2k(1 — M*)"*” / dy {sin g (4.5b) 
1], Schwarz, Berechnung der Druckverteilung einer harmonisch sich verformenden Tragflache in ebener 
Stroemung, Luftfahrtforschung 17, 379-386 (June, 1940). 
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a + E + In (4 k) ~ ran |a — cos ahe- cos ¢) + O(k’) 
yk 
me Rl(cm) 5 (I — iM’) / dy (cos ¢ — cos 2y)z’(— cos ¢) + OfK) (4.6a) 


ate 


Im(cn) = (I — M’*)™” [ dy (cos ¢ — cos 2y)z(— cos ¢) 


4) : (4.6b) 


er 


k(1 — M*)*” | dg (1 — cos ¢) sin’ ¢ z’(— cos ¢) + O(k’) 


| — 


‘ 


t 
























































F(M) = M? + In [2M7-'(1 — M’)] — (1 — M’)'? In [M™* + M“'(1 — M”’)'”] (4.72) 
5) «4 | (1 »* i) in (2M~) + 5 - (12)ae | + 0(M’) (4.7b) 
2° L\ 4° ' 2 16/° ; 
= In [2e(1 — M’)] — (3)a - yf) 
ows . 2 - 
(4.7¢) 
~ AS (1 — M’)? + O[(1 — AMf?)*] 
1) 12 
ed 025 
P| ~ 
>) | | h 
2) | | I\ 
020 t x 
| \ 
3) F(M) | | . | | 
¥ | \ 
OVS | | \ 
| \| 
) | | \ 
ta | | | | | | 
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M 
Fia. 2. F(M) vs. M, as given by Eq. (4.7). 
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where e is the base of the natural logarithms. F(.1/) is plotted in Fig. 2; for 1J > 0.85 
it is negative (having the value —0.26 at J/ 0.9), as indicated by Eq. (4.7). 

In interpreting the results (4.5) and (4.6), we observe that the real parts are simply 
the Munk, thin airfoil formulae multiplied by the Prandtl-Glauert correction factor for 
compressibility. This same statement holds for those portions of the imaginary parts 
associated with z(2), but the contribution of z’(x) to the imaginary parts is entirely 
neglected by the use of the Munk formulae; moreover the compressibility correction for 
the latter term is considerably stronger than the Prandtl-Glauert correction. 











ie J 
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U | ~ - 
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Fic. 3. Airfoil oscillating about a forward axis. 


5. Pitching airfoil. As an illustrative application of the foregoing results, we consider 
the rotation of an airfoil about an axis (P) located a distance a forward of its quarter- 


chord, as shown in Fig. 3. The displacement for unit amplitude is given by 


a | b= x (9.1) 
: 9g I h T 


Since steady-flow theory yields the real parts of the force coefficients to the desired 


accuracy, we consider only the imaginary parts. Eqs. (4.5b) and (4.6b) yield 


( a) 
Im(cia) = 2ek(1 — M?)-”%41 4 
rk(1 \! + (¢) 
(3.2) 
+(1— M°) [3 +yt+ in (3 &) - PM) |; 
= . a\* 1 h 
b\, _t ae a 
+ G- 2) ist y + wigs) — Fe) (5.3) 
a) 2 2 
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If we calculate Jm(c,,) on the basis of stationary airfoil theory, the result is 


\2 ke R ; } 2 
Im(c..) = —xi(‘) (1 — M?) vl + 1(b) | (5.4) 


For the tail on a typical airplane (b/a) will be of the order of 0.1, and k generally will 
fall between 0.1 and 0.01. Taking these values in Eqs. (5.3) and (5.4), at M = 0 and 
M 0.707, we obtain 


\ , 


ke | | 

' 0.10 | 0.01 | 

M | | | 

a _ 
Uimlemahs- | 0.83 | 0.62 

[Tentca) lises mel 


0.7 0.66 0.34 

It is evident that the use of stationary airfoil theory to predict total damping in 
pitch (horizontal surface) or yaw (vertical surface) may be dangerously non-conservative. 
Moreover, we may expect similar discrepancies in the calculation of damping in pitch 
for a swept wing. While reliable experimental evidence is small, it tends to support 
In particular, it appears that many high speed, jet aircraft exhibit a 
tendency to “snake’’, and this may well be the result of an overly optimistic calculation 
of the expected damping in yaw from the vertical fin. However, no quantitative answer 
ean be given until the questions of induction and downwash interference, (i.e. effect 
of wing on tail) are resolved. A brief indication as to the order of magnitude of the 
former effects is given in the next section. An answer to the latter question depends on 
the determination of the flow in the wake of the oscillating wing and is, therefore, rather 


this eonclusion 


more complex. An investigation of this matter is being undertaken." 
6. Three dimensional effects (Added in proof). An investigation of an elliptic wing, 
based on Reissner’s analysis of a three dimensional oscillating wing in incompressible 


and an extension of the Prandtl-Glauert rule’’, yields the result. 


flow' 
| (a\'(ac \ f /b 
Im(C, : : 9 aw _ \! a 2/5) 
h\ 9 1 3 ae 
— | (1 — M*)-'(0.530 + In A)} 1 + ( i 
( é 
ih . . a 7 \ 
= t ya — M’) | 0.01 + 0.528( ) |} 
a \ é ) 
— 0.820.M7(1 — M*)"' (6.1) 


J. W. Miles, On downwash corrections in unsteady flow, J. Aero. Sci. 17, 253 (1950). 

2). Reissner, Effect of finite span on the airload distributions for oscillating wings, Parts I [NACA 
T.N. 1194 (1947)] and IT [NACA T.N. 1195 (1947)]. 

13J. W. Miles, On the compressibility correction for subsonic unsteady flow, J. Aero. Sci. 17, 181-182 


1950 
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oo 
on 
i ¢) 


A =(1 — M*)'’@ (6.2) 


where the representative chord shown in Fig. 3 is taken as the mean aerodynamic chord, 
(dC,/da) is the three dimensional lift curve slope, corrected for compressibility, and @ 
is the aspect ratio. 

If the results of the foregoing section are examined, it is found that the dominant 
term was In k, and now is effectively replaced by —In A. It appears, therefore, that the 
two dimensional results are directly applicable only to very large aspect ratio wings, 
although the effects that we have discussed may be of some considerable importance in 
many applications. We remark that, in view of Reissner’s results for incompressible 
flow, we should expect the discrepancies between the two and three dimensional results 
to be less pronounced for larger values of the reduced frequency, so that it should not 
be inferred from our results that the two-dimensional work of Possio is necessarily in- 
applicable in typical flutter problems. 


RETARDED POTENTIALS OF SUPERSONIC FLOW* 
By J. C. MARTIN (National Advisory; Committee for Aeronautics) 


Introduction. This paper deals with a space filled with ideal compressible fluid moving 
at velocities greater than sonic velocity. A formula will be derived for the velocity 
potential at any point in space when the conditions on the disturbing surfaces are given. 
The analysis is based on the linearized partial differential equation; therefore the dis- 
turbances must be small. : 

Use is made of the concept of the finite part of an infinite integral. This concept 
was first defined by Hadamard (Reference 1). It is used to overcome the difficulty that 
arises due to the elementary solution becoming infinite on the machcone emanating 
from the origin of the source. Hyperbolic vector operators and definitions. In reference (2) 
a vector operator was used in treating source and vortex distributions in the linearized 
theory of steady supersonic flow. This and other operators are defined below. 

The hyperbolic gradient operator is defined as 
| ae ae 0 


Vh= —ifs’ —+k— l 
dx | J Oy Oz Mv 
where £ is a constant, and i, j, and k are unit vectors in the z, y, z directions respectively. 


The divergence of the hyperbolic gradient operator is defined as 


; aT ao” a” 
Veh — =< ‘a2 — a2 (2) 
Ox OY Oz 
The vector n’ is defined by the equation below 
n’ = —i6’y, + jv. + ky; (3) 





*Received November 25, 1949. This work is the basis of a masters thesis submitted to the Faculty of 
N. C. State College. The author is grateful to professors J. M. Clarkson, J. W. Cell, R. C. Bullock, and 
Jack Levine for their fruitful discussions and comments. 
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where v, , v2 , and pv, are the direction cosines of the normal to the surface to which it is 


associated. 
The divergence theorem can be expressed as 


§ E-n da - [¥-B a 


where E is a vector function and da and dv are elements of area and volume respectively. 
The normal to the element of area, da, is denoted by the vector n. The divergence 
theorem may also be expressed as 


$ W-n’ — [ va-w @ (4) 
since 
E-n = W:-n’ 
and 
V-E = Vh-W 
if 


iE, + jE, + kE, = —is’W. + jW, + kW, 


The symbol f before an integral sign will denote the finite part of the integral. 
Let E be a vector and ¢ and y be scalar functiors of z, y, z. By direct expansion the 


following identities can be proved: 


Vhy-% = Vho-Vy (5) 
Vh-yE = yVh-E + E-Vhy (6) 
If 
R= | (2-9 — Py — »)* — Be - 8)" | 
then 
v 3). vp — 2 
(3 ah = Bs @) 


Theoretical development. In supersonic flow an occurrence at a point Q(é, , £) does 
not produce immediate effects at a distant point P(x, y, z) but only during two later 
intervals when its influence is passing P. The state of the field at the point P at time ¢ 
is determined by states existing at points Q at the two earlier times (Reference 3) 


_. =O, & 

iL=t Be + Bc (8) 
_,_(@@-)M_R 

et Se 2 (9) 


where M is the stream velocity divided by the sonic velocity, c, and 6” is given by 


B? = M’-1 
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of 


The value of the potential at point P, ¢(2, y, z, t), is associated with the values 
the disturbance at points Q at times 

(¢ —&M R (2 — &)M R 

jaye wt 2 ~ SSE .S 

Bc Bc Be Be 
In the following discussion the notation ¢(x, y, z, t) will be used to denote the value 
of the potential at P at time ¢. ¢(z, y, z, t) will be divided into two functions ¢,(z, y, z, t) 
and ¢.(x, y, 2, t). d(x, y, z, t) denotes the contributions to the potential at P due to 
disturbances at points Q at times ¢, , and ¢.(z, y, z, t) denotes the contributions to the 
potential at P due to disturbances at points Q at times t, . @,(&, 7, ¢,t,) and @2(&, 7, ¢, te) 
will be used to denote the values of the potential at points Q at times ¢, and t, respectively. 

Vi 4 dp ! 4 dg k Od 


>= i wx - a 3 
9 Ox J OY Oz 


is the gradient of ¢ at P obtained by differentiating @ with ¢ constant. 


. Ob . Od; Od 
Vo, = i— —+k— 
0& On O¢ 


is the gradient of the potential at points Q at times ¢, , and 


. Id» . 92 , 1 Odo 
%. =is2+jot+kS 
0 On O¢ 


j 


ctr 


is the gradient at points Q at times ¢t, . V.¢, denotes the gradient at Q with ¢ constant. 
V ,h¢, is the hyperbolic gradient at Q with ¢ constant. 
The potential must satisfy the equation 
20 , 86 , oo 106 2V do 
re 5 a ae a oe SS 0 
Ox oy” Oz Cc ot c dxdt 


The above equation may be written 


a2 917 O* 
1d¢_2V 0¢ _g (10) 


V"hé — 3 ‘oom 
? € ay C oF Ot 
Equation (10) is a special case of one of the types of equations treated in reference 1. 
The problem is to find a function, ¢, that will satisfy given boundary conditions and the 
above partial differential equation in space. 
The hyperbolic divergence with ¢ constant of the vector 


le ( R 28 I .M a¢ 
=—Vy, — — —- —']V— ot (11 
U; R *1 Bc dt, R + * B'c Ot; 
is 
= 2V o¢ 1 0@ 
V . — X wi i — P 1 - — = — a ( 12 
hU; veh: c ox dt; c 0G ae 


The details of taking the hyperbolic divergence with ¢ constant of the vector U, , are 
given in the appendix. Likewise, the hyperbolic divergence with ¢ constant of the vector 


R 1) l ‘ _M_ db» (13) 
Bc Ot. o 


1 
on = Ge on 
‘= R ” (+. R* *BCR dat, 
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s of 
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is 





2V a6, 1 dd. 

vV . = / fA -—-= = — => =~ 

h-U, V “hz cE At, c at (14) 
Since 


o(é, n, £, 2) = o1(&, 2, £4) + be(€, 0, §, te) 


then it follows that 


Viho = Vih-U, + Vh-U, (15) 
Let 
U =U,+0, 
then by equation (4) 
$ U-n' da = [ Yau a (16) 
or 
a 2, _ 2V dd 
gu n’ da = / lv hd, co aE at, 
(17) 
_ 10%, 32 - 2¥ 2% _ 12] 
ce at; + ve c 0&dt, Cc At d 


For a supersonic flow field equation (17) becomes 
§ {bv (6-2 % R 2) 
JR Bc dt, ' Bc At, 


5 Me (2% 262) wy ” 
+ 1 3R \dt, + at,/j 9 da= 0 


ol 


(18) 


The surface of integration of equation (18) will be taken as the surface of a region 
bounded by the forecone from point P, a surface given by 
g=>2z7—-e 
where ¢ is small, and an arbitrary surface inclosed in the forward machcone from point P. 
The integral over the machcone is zero due to the concept of the finite part (Reference 
4). 
Since ¢ is small then 
iL == 
The surface integral over the area 
t=2r-—e 
reduces to a time independent problem in the limit as ¢ approaches zero. Thus this 
surface integral in the limit as « approaches zero becomes 
f 
, P11 A 0/1 
lim | | — ¢=7|s) | da = 2rd(z, y, z, b) (19) 


R an’ an’ \R 


€ 


The details of the above integration may be found in reference (4). 
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From equations (18) and (19) the integral over the surface being considered is given 
by 


_ Roo, R tt \y 1 


2 2 — Vy — ( 5 2 
mo(a, y, 2, t) + [ E o~\? — Be at, * Be Ot) R 


- Me (26 261) eo 
+ 1 POR at, + at, n’ da = 0 


(20) 


where S, , is the arbitrary surface inclosed in the forward machcone from point P. 

Since S, is an arbitrary surface, equation (20) can be applied to a region such as is 
shown in figure (1). The potential o the derivative of the potential with respect to n’ 
or both may be discontinuous across che surfaces S, . (See figure 1.) 


Az 


P (x,y,2) 




















Fie. 1. | 


For most problems dealing with linearized supersonic flow ¢ and 0¢/dn’ are zero at 
a finite distance upstream. In these cases the potential at any point is given by 


ee he ie | - ( _ R ddd, 

2, 92,9" — 3 I. E vee er a 
(21) 

R 2A¢n) t «git (246 26s) |. ta 

+ 3 at,/’ Rt BeR\ at, * o,/1" @ 
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where A¢ represents the potential difference across the surface S, . Equation (21) may 


be written 
 A°e So 2 {att 4 silt} | 
O(2, 9,2, = — 5 I li an’ an’ R da (22) 


Provided that in the differentiation in the last term of equation (22) A¢, and Ag, are 
treated as being functions of ¢, and ¢t, only. Equation (21) or (22) corresponds to the 
Kirchhoff formula, the difference being that the Kirchhoff formula applies to stationary 
mediums where equation (21) or (22) applies to mediums moving at velocities greater 
than sonic velocity. 

It is possible to represent the potential at point P by other expressions. Equation (18) 
can be applied to a region bounded by the forward machcone from point P, the same 
arbitrary surface in the forward machcone as shown in figure (1), and a surface upstream 
at infinity inclosed in the forward machcone from point P. The potential in this region 
will be denoted by ¢’. The result of applying equation (18) to the region defined above is 


, 1 ( R ae, . R oi) 1 
Se ase , cms ste oe ZS Sa 
I E o —\% — Bat, * Be at) R 








(23) 
. M dg; #44) | ” — 
tT 1BER Gs ‘qe ""* 
Equations (20) and (23) can be written as 
, ey 1 @¢ 0 i(¢ o( ty 
vannon- 2) Lebo 2 [Ogee — oo 
_—_1/f fi ag’ _— a Pot) + #4) || ' 
vise 2r [ it an”’ an” | R a (25) 


where in the differentiation with respect to n’ or n”, ¢,(t,)do(te) , di(t:) and 3(t,.) are 
treated as being functions of t, and t, only. The result of adding equations (24) and 





(25) is 
ee eee et (26, ) 
9(2, ¥, 2, 0 = 2r I. th an’ + an” 
(26) 
_ a EE + da(t2) = o%(4) — sa) |} ‘ae 
an’ R 


Since ¢’ is an arbitrary solution of equation (10) then it may be chosen so that either 
0¢/5n’ + d¢'/dn” or 





x [a + do(ts) — o4(t) — #56) | 
on’ R 


is zero; therefore ¢(x, y, z, t) may be given by 


1 1 /d dg’ 
¢(z, Y, 2, t) oe R (36, + ov) da (27) 


2r / 8, 
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or 
be a th) + do(te) — b1(ti) — $3(te 

=a | _ ; | J + Gol 2) — di(ts) — $2 2] da (28) 

nr Js, ON R 

depending on whether ¢’ is chosen so that 0¢/dn’ + d¢’/dn"” or 


é k (t,) + do(ts) — $a) — ota 


an’ R 
is zero. It should be remembered that in the differentiation with respect to n’¢,(t,), 


$2(t2), di(t,) and $3(t.) are treated as being functions of t, and ¢, only. 


APPENDIX 


The evaluation of V,h-U, . The valuation of V,h-U, requires the use of equations 


(5), (6), (7), (8) and (9). From equation (6) V,h-U, can be expressed as 


l R Oo 1 1 ° M Od 
. = . _ ‘a . > —_ . : - 
v h-U, Vih R V4, Vv h-¢,5 R + Vh 3c al, v R + Vih i eR at, (A-1) 


The first term of equation (A-1) may be expanded to give 

l l M l Op l 
Vg, = ee ahs me RS ee ek 

or cR O€ At, 8°cR ' Ot; + “o R 


Vh. = 
h R I 


> V hd; — 


The second term in equation (A-1) may be expanded to give 


l l l . M do 1 O¢ 
ain Vi an «Ve. Vax 2 0) eS 
. 1 R R hg, j R : Cc Ot; Rc Ot, 


The third term in equation (A-1) can be written as 


Vi : R 0d, Vv 1 = 1 ; 0g; aie. ke Vv l ; Vh Od, om RM Vv L he od: a - od: 
Bc ot, R cn bt, Be &£ Ot; ro R Ot; BcR dt; 


The last term in equation (A-1) can be expanded to give 


Vhei M 0¢,_ _ M 4, 9 ( 1} M dq _ M 0%,  —M_ 0&9, OR 
rl BR dt, c 0t, 0¢\R cR 0é dt, BCR Ot; BcR dt; O€ 


From the above results equation (A-1) can now be written as 


l 2V a “b | O° 
TA a — Os, — =: 1 — —, —*} 
U, R ; , Re O& Ot; Rec Ot; 


which was to be shown. 
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NOTE ON THE STRESSES PRODUCED BY NUCLEI OF THERMO-ELASTIC 
STRAIN IN A SEMI-INFINITE ELASTIC SOLID* 


By BIBHUTIBHUSAN SEN (Presidency College, Calcutta) 


1. Introduction. It has been shown by Goodier' that if an element of volume dQ 
surrounding a point C in an infinite elastic solid is at a temperature T while the rest 
of the solid is at zero temperature, the displacement at any point is the gradient of a 
function 


a(1 + ») T : dQ (1) 


(l—v) 4 r’ 


where v is Poisson’s ratio, a, the coefficient of linear expansion and r, the distance of 
the point from C. The singularity at C is the simple ‘centre of dilatation’ which is supposed 
to form a nucleus of thermo-elastic strain at the point. The problem solved in this paper 
is that of determining displacements and stresses in an isotropic semi-infinite elastic 
solid which has a nucleus of thermo-elastic strain as defined above at any point C below 
the plane boundary. 

2. Statement of the problem. We suppose that the solid is bounded by the plane 
Z = 0, the axis of Z being drawn into the body. The coordinates of C are taken as 
(0, 0, c). Putting 


_al+yv) Td (2) 


>? = ‘ 
v= (1 — pv) 4r ’ 


and 


3 L Se ee 
iia RQ, [fe +y¥+e-o0')]”’ (3) 


we can express the components of displacement wu’, v’, w’ due to the nucleus of thermo- 
elastic strain situated at C as 


es yf = Pp wo’ = P > (4) 


u’ = | ar? ay ae 


These components give us the stresses 


i _EP dv ai _3EP _ x{z — ¢) (5a) 
(l+v)0xdz (1+) R° 
7p Pa 7p —— 
__ - x. v, _ _3EP_ ye — ¢) (5b) 
(1 + v) dy dz (1 + ») R? 
EP av EP | 3(z - | 
‘—- —— a = So —_— —F SS i 5 
7 a+na  (1+>) | Rit Re (Se) 


E being the Young’s modulus. 
*Received Dec. 19, 1949. 
1J. N. Goodier, On the integration of thermo-elastic equations, Phil. Mag. 23, 1017 (1937). 
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In order that the surface Z = 0 may be free from stress, we must superimpose a 
stress system 72! , r/! , a!’ possible in the body such that 
3EPxe : 
[727] - —([ri.], o> ’ (6a) 


- (1 + y(a + y? + 2) 5/2 


3EPyc 


[tre]emo = —[trele-o = — rg age (6b) 
(1 a v)(a- + y 4 r+ iad , 
EP 1 BY 
[o;’ ae = a: z=0 — 7 | 2 2 wa «= lho: |. ~ “ =| (6e 
: dia (l+y)L@ +y +c)” (fo +y +c)” 
These stresses 7!’ , 7)? , of’ can be obtained easily by the following method which was 


given by the author in a previous paper.” 
3. Solution. Since by assumption r!! , rj! , o!’ are stresses possible in an isotropic 
elastic solid, they must satisfy the equations 


1 0°06. 


Pads = 0 75 
at (1 + v) Ox dz : (7a) 
, l 0, 
27 : 2 = 7} 
Vite + (1 + v) dy dz : (75) 
I ” 0, " 
Viel +a+y a ~% (70) 
where V7’ is the Laplace operator, and 6, = a2’ + of’ + ot’. 0, being a harmonic func- 
tion, we can obtain the solutions of the above equations as 
ia I 00, 
ee "Si +7 * Ox T ot» ™” 
” I I. : 
Te = —50 +) * ay + ds, (8b) 
1 30, 
o;' = 30 +» 2 + ¢3, (8c) 


where ¢, , ¢2 and ¢; are harmonic functions. Using the results (6) and (8) we find that 


_ _BEPxc (9a) 


[¢:].-0 id [es z=0 = (1 + v)(20? ry y? 4+ ¢)”? ’ 
2D Dew 
B3EPyc __ (9b) 


daleeo = [Te]-0 = ; ' 57/2? 
Ie] Q+r(a@ +y +c)” 


EP l 3c” 
[3].-0 = CAg a so re os 5 Oe So. a ee (9c) 
(l+y»yL@ +y + ¢) (2? + y? +c’) 
2B, Sen, Stresses due to forces and couples acting in the interior of a semi-infinite elastic solid, Bull. 
Cal. Math. Soc. 32, 73-83 (1940). 
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Putting 


I l 
>o=—- ; > Fre. a ee , f. § _reaueee 10 
me R, [2 +y + (e+0)']"” (10) 


we find that the conditions (9) are satisfied if we write 





EP 0°vp EP 0°v. EP av 
~i = = “2 d: = 7a —-2, :=-7, ——§. (11) 
(1 + ») dx & (1 + ») dy & (1 +») dz 
Since 
OT: O74: dc,’ a 
Ox Oy Oz 


and v, is a harmonic function, we have from (8) 


30, 8d, , db» os | wn OD 
_- =| = —4EP -> 
Oz ai + n| 2 * Oy + Oz nd dz*’ 
whence we obtain by integration 


2 
OV. 


aa? (12) 


6, = —4EP 


The function of integration in this case as well as in the cases that will follow, will be 
zero since the stresses and components of displacement must tend to zero as Z tends to 
infinity. 

Substituting into (8), the values of ¢, , @2 , ¢; and 6, obtained in (11) and (12) re- 


spectively, we have 


IP 93), 20, 
i | 20 5 + Z| (13a) 


(1 +») L™ dx dz” * ax dz 
CE ES) eT (13b) 
a = Ea me - FI A309 
Combining these results with those given in (5), we get 
ra tat oth = GA | B+ mS (14s) 
Tye = Te + Te = ( . ES +S eit + Zu | (146) 
s, wet eee rer | + 22 oe - a. (140) 


If of , o, , of be calculated from wu’, v’, w’ given in (4), it will be found that 


6i,=¢c+04+05=0. 
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Hence 


OV». 


02° * 


A=: A. + é. = A, = = {EP 





Now let wu, v, w be the components of the resultant displacement parallel to the 
axes of x, y and z due to the two systems of stresses with one and two dashes respectively. 





Since 
Ow ee Ov OV. a°v. Ov 
_ | ) ) 1 bt oe = 3 
Oy Lies + vho, — v0) = PL Oe 4 a Se Hy gg MH 
Oz EK s ; Loz Oz Oz ad Oz 
we find by integration 
ov . oO OV> 
w= P| — + 22—= -— (3 —-— 4) — 
oz Oz Oz 
(16) 
aT dQ(1+p [ —C fyo(z + ce) — (2 + 3c) Gz(z + c)” 
dor (l — p | R ig R i 
f 
A ' 
Again 
Ou 2(1 + v) ow 
dz E  @ 
Ov; Oo —— OV. Ol ] 
= P| — 22 — +. i? 2. 4. (3 ty) 
OX Oz Ox O. Ox Oz On dz _| 
By integration we get 
Ov, ra) rel) 
“ue Fi — + 2 + (3-4 
Ox Ox Oz C2 
(17) 
aT dQ (1 + | 1  (3—4) Gele+ 2 | 
= ? a ; - e ~ ‘ 
dor (1 — LE is R 
Similarly from the relation 
On 2(1 + pv) Ow 
dz E Oy’ 
we have 
ef dQ(i+v) [1 , (@—4y)  6ee+ | ae 
oa SS Fe Yl p3 . 3 SS . 3) 
Ar (1 — pv)‘ LR’ ie i. 


Components of displacement being given by (16), (17) and (18), the stress com- 
ponents (other than 7,, , T,. , @: which have already been obtained in (14)) can be 
easily calculated. 

We find that on the plane Z = 0. 
aT dQ (1 + ve 


[w].-0 = = Re (19) 





e 
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where 
Re=ae?t+yte’. (20) 


It is evident that if there be a continuous distribution of such nuclei (as enunciated 
above) over a certain region, we can find the total effect by integration over the region. 


PRACTICAL SOLUTION OF SIMULTANEOUS LINEAR EQUATIONS* 
By O. L. BOWIE (Watertown Arsenal, Watertown, Mass.) 


I. Introduction. The effect of “rounding off” errors in the solution of the set of 
simultaneous linear equations with constant coefficients 


n 
Saye, =e = 1,% 0050) (1) 
= 

can lead to considerable inaccuracy of solution. Furthermore, the accuracy of an approxi- 

mate solution of (1) usually cannot be determined from the errors in the individual 

equations without an excessive amount of additional computation. Unfortunately, most 
conventional techniques lead unavoidably to “rounding off” errors and the computer is 
often faced with both of the above problems. 

When the form of (1) is such that the well-known’ “Classical Iterative” technique can 
be applied, the difficulties above are usually not encountered. In most practical problems 
the coefficients and the desired solution can be handled well within the digital capacity of 
the average desk calculator. Since only the coefficients and approximations of the solution 
occur in the process, the necessity of “rounding off’’ is therefore usually eliminated. Accu- 
rate estimation of the accuracy of an approximate solution is inherent in the iterative 
character of the process. 

Although rigorous investigations have been made as to the form of (1) sufficient to 
ensure convergence with the Iterative technique (e.g. footnote 1), the criteria obtained 
are too impractical to be applied for each case. For all practical purposes, however, a 
simple rule of thumb is usually sufficient to consider. If (1) can be so ordered that 


Q;; > (aj, 9 Bye 9 °** 9 By sat » Big4a 5» °** » jn) 7 - 1, 2, "o® 2 Oe (2) 


it can be safely concluded, in general, that the Iterative technique will converge. 

It is the aim of this note to demonstrate how conventional techniques can be prac- 
ticably utilized to transform (1) into an equivalent system which satisfies the conditions 
(2), thereby enabling the computer to take full advantage of the digital type of solution 
possible by the Iterative technique. It is the opinion of the author that in general, it is 
more practical computationally to regard techniques involving “rounding off” errors such 
as analogue computers, Crout’s Method,’ ete., as auxiliary tools for obtaining an equiva- 
lent system satisfying (2) rather than methods for obtaining a final solution directly. 


*Received Dec. 30, 1949. 

1H. Hotelling, Some new methods of matrix calculation, Ann. Math. Stat. 14, 1-34 (1943). 

2P. D. Crout, A short method for evaluating determinants and solving systems of linear equations with 
real or complex coefficients, Presented at A.I.E.E. Summer Convention, 1941. Published by A.I.E. 1941. 
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II. A procedure for the solution of simultaneous linear equations. Let us consider 
systems of equations (1) which in their original form are such that the Iterative technique 
leads to a process which diverges or converges too slowly to be practical. Systems of 
equations (1) which have unique, non-trivial solutions can be represented in matrix form 
as 

AX =(C, (3) 


where A is the n X n matrix (a;;) of rank n, X is the n X 1 column matrix (z;), and C is 
the nm X 1 column matrix (c;) where C # 0. Although a given system of equations (3) 
will have a unique solution X, the converse is not true. In fact, let Q be any n X n matrix 
of rank n. Then the system of equations 


QAX = QC (4) 


has the same solution as (3) or is “equivalent.” 

*‘Diagonalizing”’ operator matrices, Q, can be considered as those operators for which 
QA consists of non-zero elements on the principal diagonal alone. Obviously the inverse 
matrix, A~', belongs to this class since A~'A = J, where J is the identity matrix. Little 
would be gained computationally if, for the solution of a given system, the exact calcula- 
tion of a diagonalizing operator were required. Such a calculation would be subject to 
the difficulties arising from “rounding off”’ errors. 

On the other hand, the notion of “equivalent”? systems when combined with the 
Iterative technique lends itself to a practical procedure for the solution of systems of 
simultaneous linear equations. An approximate diagonalizing operator will in general, 
be sufficient to determine an equivalent system of such a form that the Iterative technique 
can be successfully applied. ‘‘Pseudo” diagonalizing operators, Q, will henceforth be 
considered as those operator matrices such that QAX = QC is in a form such that (2) is 


satisfied. 

Algebraically, the elements q,;; of a pseudo diagonalizing operator matrix can be 
determined as follows: Let D = (d,;;) be ann X n matrix of rank n such that the elements 
d;; satisfy the type of condition (2). The elements g;; must be determined so that 


(9:3) *(a;;) - (d;;) (i, j = l, 2, etic , 2). (9) 


Equation (5) yields the n systems of equations for the successive rows of Q, 


» Aides = dys (7 = 1,2, °-: , 7) (6) 
j=1 
corresponding to k = 1, 2, «++ , n. 


For example, if d,; is chosen as the identity matrix, then the systems (6) represent the 
algebraic conditions for the elements of the inverse matrix, a 

III. Practical calculation of pseudo diagonalizing operators. It should be noted first 
of all that the solution of (6) by most conventional techniques requires little more work 
than for a single system. This is true because the coefficients in (6) are the same for each 
of the n systems. 

In analogue computers, the successive rows of the matrix q,; are economically obtained 
by setting up the coefficients a;; once only and then adjusting successively the right hand 


sides in the order indicated in (6). A practical procedure for the solution of (6) on auto- 








No. 4 


sider 
ique 
is of 
‘orm 


(4) 


ich 
Tse 
ttle 
ila- 

to 


the 
of 
‘al, 
jue 
be 


is 


d 
d 


)- 








1951] 0. L. BOWIE 371 


matic calculating devices is obtainable by an obvious extension of “Crout’s Method.” 
Since the coefficients in (6) are the same for each system, an obvious extension of ‘“Crout’s 
Method” would be to carry along (n — 1) additional columns corresponding to the vari- 
ous right hand sides of (6) in his “auxiliary matrix.” 

The most practical consideration in the solution of (6) is the degree of freedom the 
computer has in the choice of d;; , the only restriction being that d;; satisfy the type of 
condition (2). In analogue devices, the computer could adjust iteratively the values of 
q:; for example, until a suitable set of d;;’s is obtained. Since the original draft of this 
note was written, this possibility was investigated’ and led to very successful results. On 
the other hand, if techniques designed for the automatic calculator are employed, the 
most logical approach would seem to choose a particular d;; , (for example, the identity 
matrix) and then examine approximate solutions of (6) so obtained to see if a suitable set 
of d;;’s occur correspondingly. 

Since the dominance of the diagonal coefficients indicated in (2) need not be too great 
in practice, in fact ratios of 10 to 1 are often sufficient to ensure rapid convergence of the 
Iterative technique, it is evident that quite crude solutions of (6) will often yield pseudo 
diagonalizing operators. In fact, if the system (1) is reasonably “well-behaved’’, one can 
fruitfully consider solutions of approximations of (6). In the typical problem carried out 
in this note, the coefficients a,; were drastically rounded off in the determination of pseudo 
diagonalizing operators, yet reasonable values of q,; are obtained. 

If the identity matrix, 7, be chosen for D, as was suggested for application of tech- 
niques designed for the desk calculator, one might consider a pseudo diagonalizing opera- 
tor as an approximation of the inverse matrix, A~'. It is interesting to note, however, that 
whereas a good approximation of A~' can be usually considered a satisfactory operator, 
the converse statement is often far from true. 

When several of the equations of a given system are of a form satisfying the condi- 
tions (2), calculation can be shortened by carrying these particular equations into the 
equivalent system directly. To carry the p’th equation of the original system unaltered in 
form into the equivalent system, the p’th row of the pseudo diagonalizing operator can 
be immediately written as zero for all elements except that corresponding to the p’th 
column, and the element g,, is written as unity. 

IV. Pseudo diagonalizing operators from inspection. A technique will now be out- 
lined by which operator matrices for improving the rate of convergence can be determined 
by inspection. 

Consider the systems of equations 


RX = 5S, (7) 
where the matrix R can be expressed as the sum of a diagonal matrix of the form d/J and a 
remainder matrix, ¢, such that lim ¢{... = 0. In the matrix product d/, I is the identity 
matrix and d is a scalar greater than or equal to unity. Thus, (7) can be written as 
(I+e0o)X =S. (8) 
Operator matrices which yield equivalent systems for which the rate of convergence is 
improved can be written by inspection. As examples, 


88, Ackerman, A low-cost simultaneous equation solver, Watertown Arsenal Report Number WAL 
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Q, = (dJ — © = (2dI — R) 
0d, = (d’] — de+ €) (9) 


etc. 
yield the equivalent systems, respectively, 


(FI — &)X = (aI — 08 


(d°I + e)X = (I —de+)S (10) 
etc. 
Since lim &§_.. = 0 is part of the hypothesis, it is evident that, in general, the domina- 


tion of the diagonal coefficients in (10) will be increased as compared with (8), thus leading 
to more rapid convergence in the solution of the equivalent system. The operator Q, would 
seem the most practical for it can often be written down mentally. 

V. An Example. The example which is chosen to illustrate the procedures in ITI and 
IV is one of the systems used to illustrate “Crout’s Method” in footnote 2. Consider the 


system 


Zy Ls Is Le = C¢ 
12.1719 27.3941 1.9827 7.3757 = 6.6355 
8.1163 23.3385 9.8397 1.9474 = 6.1304 (11) 
3.0706 13.543: 15.5973 7.5172 = 4.6921 
3.0581 3.1510 6.9841 13.1984 = 2.5393 
The calculation of a pseudo diagonalizing operator Q was made from (6) with D = | 


(identity matrix). A ten bank Marchant calculator was used in conjunction with the 
extended “Crout’s Method.” The coefficients a,; were rounded off to one decimal place 
and a maximum of three significant figures were kept in the subsequent computations. 


The pseudo diagonalizing operator thus obtained is 


(—1.50 2.80 —2.00 0.91| 
0.56 —1.00 0.75 —0.37| 

Q = |-0.44 0.76 —0.46 0.23 | (12) 
0.44 —0.81 0.53 —0.17) 


The corresponding equivalent system, found by multiplying (11) by Q becomes 


Zi Z2 Zs nr" 
1.109461 .037260 — .261959 — .234686 = .138433 
— .128583 .993876 . 384470 — .062516 = .165014 
.103639 .178622 1.037369 .092436 = .165157 (13) 
— .111026 — .208449 — .018497 .978302 =  .009128 
”? 


The system (13) was solved by the “Classical Iterative” method to nine decimal place 
accuracy in ten iterations, yielding 


(.15929 1130,  .14691 7740, 11257 4804,  .06084 0731) (14) 
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It is interesting to compare (14) with the solution obtained by Crout. In footnote 2, five 
significant figures were kept in the computation and the solution of (11) was given as 


(0.15942, 0.14687, 0.11261, —-0..060806.). (15) 


The equivalent matrices would be in better form of course, if a more careful calcula- 
tion of (12) were made. Thus, keeping five significant figures in the solution of (6), a 
pseudo diagonalizing operator was obtained for which the equivalent system was solved 
to nine decimal accuracy in four iterations. 

To illustrate the methods of Section IV, an approximate value of Q, in (9) was chosen 
to accelerate the solution of (13). To maintain a reasonable number of significant figures, 
the coefficients of (13) were mentally rounded off to two decimal points. By inspection, 


therefore, the improving operator, Q, , becomes 


( 39 —.04 26 23) 

13 1.01 — U6 | 
Q,=|-—.10 —.18 96 —.09 (16) 

ll 21 02 1.02) 


The equivalent system formed by multiplying (13) by (16) required only five iterations 
to yield (14 

ACKNOWLEDGMENT: The author is indebted to Lyman B. Stewart, Watertown Arsenal, 
for examination of the computational feasibility of the procedures. 


PLASTIC DESIGN OF BEAMS AND PLANE FRAMES FOR 
MINIMUM MATERIAL CONSUMPTION* 


By JACQUES HEYMAN (Brown University)** 


Abstract. This paper is concerned with the design of plane frames in such a way 
that the material consumption is a minimum. The method of solution is to set up linear 
inequalities for the variables involved, and to solve these inequalities by the Dines' 
method. Three slightly different classes of problem are treated; collapse design under 
fixed loads, collapse design under varying loads, and shakedown design under varying 
loads. Illustrative examples of each are given. 

1. Introduction. The problem discussed in this paper may be stated as follows: 
Given the geometrical configuration of a frame which is acted upon by specified loads, 
how should the various members of the frame be proportioned so that the total ma- 
terial consumption is as small as possible? The type of frame to be considered consists 
of straight members, all in the same plane. Between joints, the assumption is made 
that any one member has a uniform cross-section. The moment-curvature curve for any 
cross-section will be taken to have the general form shown in Fig. 1, where a limit 


*Received March 7, 1950. This paper is a shortened version of Technical Report A1l-45 prepared 
for the Office of Naval Research under contract: N7onr-358. 

**Now Fellow of Peterhouse College, Cambridge, England. 

1. Dines, L. L. Systems of linear inequalities. Ann. Math. (2) 20, 191 (1918-1919). 
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moment M is attained at a finite curvature such that any further bending is accom- 
panied by no increase in M. This limit value M is known as the full plastic moment, 
and the section at which it occurs is called a plastic hinge. It will be assumed that the 
joints connecting the members of the frame are sufficiently rigid for the full plastic 
moments to be developed (if necessary) at those joints. 

The total material consumption x for a given frame is given by the proportionality 


za > M,l,, (1) 


1=1 


where M, is the full plastic moment of the 7‘* member of the frame, and J; is its length. 
This linearization of the material consumption is not strictly correct, but for the pur- 
pose of developing suitable methods for design, it will be assumed that a continuous 


linear range of sections is available. 


MOMENT 











CURVATURE 


Fic. 1 


The general state of a frame of n redundancies can be expressed as the sum of one 
arbitrary fixed equilibrium state and n linearly independent residual states. By a state 
is meant some distribution of bending moment, so that a state-in equilibrium with the 
applied loads is any bending moment distribution such that equilibrium is attained. 
A residual state satisfies equilibrium when the applied loads on the frame are zero- 
Thus, confining attention to any one cross-section in a frame, the bending moment 
there may be expressed as 


M*+e¢Mi+eMi+---+c¢,M, (2) 
where M* is the equilibrium bending moment at the section, and M{ M? are the 
linearly independent residual moments. Suppose the full plastic moment at the section 
(as yet undetermined) is M. Then since the total bending moment cannot exceed the 
value M, 

—-M< M*¥+e¢Mi+eoMi4+---+¢M, <M (3) 


and continued inequalities similar to (3) may be written for every critical section of 
the frame. Only concentrated loads are considered in this paper, so that shear forces 
are discontinuous at the joints and under the loads, and are constant elsewhere. This 
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means that the bending moments in the frame are linear, and critical sections can occur 
only at the joints and under the loads. Writing down the inequalities (3) for these 
sections, the problem may be solved by Dines’ method, and the ideas are best presented 
by means of examples. On account of the limitations of space, these examples will deal 
with the design of beams; the design of plane frames introduces no new ideas but in 
general entails rather more arithmetical labour. (In the full report, of which the present 
paper is a shortened version, a rectangular portal frame is treated). 

2. Example 1 (fixed loads). Consider the two-span continuous beam on three supports 
carrying the central fixed loads P; and P, (P; > P:) as shown in Fig. 2. Since the sys- 


a 


Fetatatel 


Fia. 2 





tem has one redundancy, the plastic behaviour can be represented as the sum of an 
arbitrary fixed equilibrium state and one residual state. These may be taken as the 
two bending moment distributions in Fig. 3, where the factor c corresponds to the 











: ' C 


(b) 


Fic. 3 


arbitrary multiplying factors on the residual moments in expression (2). Denoting the 
full plastic moments of the two spans by M, and M, (the values being as yet undeter- 
mined), the condition that the total bending moments at the critical points should be 
less than or equal to the full plastic moments at those points gives the following in- 


equalities: 


(1) under the load P, , —-M,<pt+te<sM,, 
(2) at the central support, —-M,<2<M,, 
(4) 
—-M,<2<M.,, 
(3) under the load P, , —-M,<pm+e<sM,, 
where p, for example is given by p, = P,l/4, the “free’’ bending moment shown in Fig. 


3(a). The two inequalities for the central support are obtained by considering respec- 
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tively points slightly to the left and to the right of the support, since it is not known @ 
priore whethe1 M, S M, ‘ 


For the purposes of calculation, each continued inequality in (4) may be expressed 


as two simple inequalities, and the set (4) rewritten: 


c+p+M,>0, -c-—p+M,2>0, 


M M, 


+2120, -c¢ +720, 

A M, 
+*2>0, -c +720, 
c+p.+ M, > 0, —c—-p,+M, > 0. 


For the present purposes, the actual value of ¢ in the set (5) is of no interest, so long as 
it is known that a value of ¢ exists such that each inequality can be satisfied. Dines 
has shown that necessary and sufficient conditions for the existence of such a value of 
c are obtained if every inequality in the set (5) with a coefficient of +1 for ¢ is added to 
every inequality with a coefficient of —1 for c. If these additions are performed, c¢ will 
be eliminated, and the resulting set of inequalities (4 * 4 = 16 in this example) pro- 
vide the necessary and sufficient conditions for the existence of c. 

In the actual operation, it is found that some inequalities are redundant; for ex- 


ample, adding the first and second of set (5) gives 17, > 0, while the second and third 
give VW, > 3p, , so that IJ, > O is redundant. The five remaining inequalities out of 


the possible sixteen are 


— Pr + 2M, > 0, 

— Pp + M,+ 4M, > 0, 

—-p+p+ M+ M,>0, (6) 
—p,+4M,+ M,>0, 
— p + 3M, > 0. 


The material consumption parameter x will now be introduced into set (6). Since the 


two beams have equal spans, 7 may be taken as (see 1)) 
zs=M,+M.. (7) 
M, may be replaced in set (6) by (2 — J/,), giving, upon slight rearrangement, 
— M.+ x-— 3p, > 0, 


— M, + 2x — 2p, > 0, 


M.+ x — 2p, > 0, 


La) 

ee 
t 
wold 


— 4, > 0, 
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together with 
x > (~pi — Pr). (9) 


Now for the problem of determining the minimum value of z, the value of M, is not 
required, and Dines’ method may be employed again on set (8) to eliminate M, , giving 


z> tpt P2 


ra 3P1 + 3P2 » 
(10) 


8 
IV 


3): + 3p2 » 


zr> p+ 3p. 
Inequality (9) and the second of (10) are obviously redundant, and for p, > p,. (as 
assumed), set (10) reduces to 

t>pit ip, (11) 
the other two inequalities also being redundant. This single inequality (11) is a necessary 
and completely sufficient condition that values of 1, , /, , and c can be found to satisfy 
the original set (4). Since it is required that x should be as small as possible, the equality 
sign will be taken in (11), so that 

z= pi + nr (12) 


that is, the existence condition is just satisfied. It is to be expected that unique values 
of M, , M. , and ¢ will be derived, and this is in fact the case. Inequality (11) was ob- 
tained by the addition of the second and fourth of set (8); the second of (8) states that 
M, < 2p, , for the minimum value of x given by Eq. (12), while the fourth gives M, > 


3p. . Hence 
M, = $2; M, = pi — 4pr2 (13) 
and, by going back to the original set (5), 
2c = —3p, = —M.. (14) 


The bending moment distribution resulting from the analysis is shown in Fig. 4, which 


[Bx Ke73 
NY 





M,=(P,- SP) 
Fia. 4 
shows that plastic hinges are formed at all three of the critical points. Physically, this 


number of hinges is sufficient for the structure to become a mechanism and hence to 


collapse. 
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3. Example 2 (variable loads, collapse design). Consider the same beam shown in 
Fig. 2 but with the loads varying between the limits 


1<P<Q, —-Q% <P,2<5Q, 12>. (15) 
The same set of inequalities (5) may be written down containing the parameters p, 
and p, . The sixteen inequalities resulting from the elimination of ¢ must be examined 


in the light of the fact that p, and p, vary, and the worst values of p, and py (i.e., +4q: , 
tg.) inserted in each inequality. After removal of redundant inequalities, the set 


becomes 


= a + aM, > 0, 
Oh + M,+ 4M, > 0, 
— 0; go+ M,+ M,2>0, (16) 


Operating on this set as before to find the minimum value of z, it is found that the 
third of (16) is the controlling inequality, giving 

L=Qat aq. (17) 

This is one of the few examples where unique values of M, and J/, are not generated 

by working back through the successive sets of inequalities; instead, ranges are found 


Aas follow Ss: 


(18) 
2q. > M, > 4a, 
(3% +g) 2 M,. 
As a specific example, suppose gq, = g2 = g. Then 
} 2 4 2 
M, + M, = 2g, 3 9 ~My 3 4) 3 4 > Ms > , 4% (19 


and any values of 1, and M, satisfying (19) will give a constant material consumption. 


It is perhaps of interest to note that if 2 = E(.W7"l), wheren < 1, the minimum material 
consumption is given by M, = 2M, = 4q (or vice versa), the worst case occurring for 
M, VU q. An asymmetrical solution is obtained for what appears to be a com- 
pletely symmetrical problem. In practice, for n = 0.6, which is a reasonable value for 


the exponent, the worst error due to the linearized material consumption amounts to 
less than 2% in this example!. 

4, Example 3 (variable loads, shakedown design). If, after a certain arbitrary time, 
a structure subjected to varying loads resists all further changes of load in a purely 


in 


Pi 
ed 


6) 
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elastic manner, then the structure is said to shake down under the given loading system. 
It may be shown* that a necessary and sufficient condition for shakedown to occur is 
that a sum of residual states can be found such that when it is added to the elastic 
response of the structure, the moment at any cross-section is less in magnitude than the 
full plastic moment at that section, that is (ef. Eq. (3)), 


—-M<M,+¢,M{ +c.M;, + ---c,Mi < M, (20) 


where J/, is the elastic solution for the problem under the given loads. Now, as the loads 
on the frame vary, M/, will vary between limits, say M?* and M?"*, so that (20) may 
be replaced by the two inequalities 


Me™ + ¢,Mi + o.M, + --- +¢,Mi < M, 
(21) 
— M < MP" +¢,M/+c¢Mi+---+cM. 


The elastic solution of a given structure presents no formal problem, although for 
highly redundant frames the numerical labour may be heavy. The bending moment at 
any cross-section may be expressed in terms of the given loads and the flexural rigidity 
(B) of the various portions of the frame, so that the elastic parameters (flexural rigidities) 
used in the solution are not the same as the plastic parameters (full plastic moments). 
However, once a certain class of cross-section has been decided upon, an empirical 
relationship may be plotted connecting the two, so that if a full plastic moment is given, 
it is simple to read off the corresponding flexural rigidity of the section, and vice versa. 

In order to work the following example analytically (and it should be stressed that 
any given relationship will not lead to more complicated results), it will be assumed 
that 


B= kM*”, (22) 


where / is a numerical constant. The factor of 4/3 expresses the fact that approximately 
geometrically similar sections are involved. 

Consider the same beam as in the previous examples (see Fig. 2) with the loads 
P, and P, varying arbitrarily between the limits 


os Ff. SF; 0< P, < 2P. (23) 


If the flexural rigidities of the two spans are B, and B, , the elastic bending moments 
at the critical points are: 


5B, 8B,) — pl: 
(1) under P, , Ma, a5 + 8B,) - pA3B,) | 


8 (B, + B) 


3B, (3B, 
(2) at the central support, M,: = ~1f 2088) + GB) | (24) 
(3B, (5B 8B. 
(3) under P, , M,; = i - Pil Ba) + m(SBy + BB) | 
« \¥71 2 


2Neal, B. G., The behaviour of continuous beams and plane frames under repeated loading, Technical 
Report Al1-32 to the Office of Naval Research under contract N7onr-358. 
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Hence as P, and P, vary between their limits, 


vo — of $8: + 85, | | ~ — ~# 9 | 
~ 7* “SL B, +B. ) ve SLB, + B,J)’ 

7 ab, + 98, - 
Ma” = 0, Ma = of B,+B, 1’ ated 
aie of 108, + 6B | ae 38, | 
ee ~ = ee FF ~— 48.4 8,7 


These values of moments are non-linear in the values of B even for this simple example, 
and in any case the relationship between B and M (the full plastic moment) is not 
expressible analytically in general, so that an iterative numerical method will be de- 
veloped. The general method of solution may be tabulated as follows. 
1. Assign flexural rigidities to the various members of the frame, specifying these 
rigidities in terms of one unknown. 

2. The elastic solution may be obtained, and the limits of the bending moments at 

critical sections calculated as the loads on the frame vary. 

Forgetting that flexural rigidities have already been specified, use these elastic 

limit moments as a state of the frame to which must be added some sum of residual 

states such that the bending moment at any section is less in magnitude than the 
full plastic moment. 

4. Set up the usual inequalities, and solve the problem by Dines’ method, which 
will give the actual full plastic moments necessary for minimum material con- 
sumption. 

5. From these full plastic moments, calculate (either by means of an empirical curve 

or by some relationship such as Eq. (22)) the flexural rigidities of the various 

members. If these rigidities are far removed from those originally assumed, 
obtain an elastic solution with the new values, and rework the problem. 
Continuing with the two-span beam example, if it is assumed that 
B, = 2B, (26) 

(since the right hand span has to support twice as much load as the left hand), then, 

working to the rules given above, it is found that 


t= i p = 2.33p, M, = =p = 0.67p, M,= 2 p = 1.67p. (27) 
Hence for the next trial, from Eq. (22), take 
; = (He) = (2.5) = 3.4. (28) 
This ratio of flexural rigidities gives 
z= 2.68p, M,=0.78p, M, = 1.90p, (29) 


and the next trial, although B,/B, , is reduced to 3.25, gives these same values of M, 


and M,. 
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5. Conclusion. The methods presented above have been applied to more complicated 
frames and readily give the required solutions for collapse design under constant or 
varying loads. For shakedown design, the iterative numerical method converges fairly 
rapidly. It may be found easier for highly redundant frames to obtain a new elastic 
solution at each stage using the numerical values obtained from the previous analysis, 
since an analysis with numerically unspecified flexural rigidities is extremely tedious. 
lor all three types of design, the introduction at an early stage of the numerical values 
of the loads simplifies the work greatly, since it will be found that a large proportion 
of the inequalities generated become redundant and can be ignored. 

Only examples of concentrated loads on straight members of uniform cross-section 
between joints have been examined, making it possible to pick by inspection the critical 
cross-sections. However, the basic ideas are not altered by the introduction of other 
variables; the analysis will be more complicated, but aids to calculation may be in- 
troduced which leave the basic problem unchanged. 


THE METHOD OF CHARACTERISTICS APPLIED TO PROBLEMS OF 
STEADY MOTION IN PLANE PLASTIC STRESS* 


By P. G. HODGE, JR. (University of California at Los Angeles) 


A method is outlined for obtaining the stress, strain, and thickness distribution in a 
thin sheet which is strained plastically in its plane. For the particular case of steady 
motion, a method is given for obtaining directly the final solution to certain types of 
boundary value problems. A step by step procedure is indicated for the general case of 
non-steady motion. 

1. Introduction. This paper is concerned with the stress and strain distribution in 
a thin sheet which is strained plastically in its plane, under conditions of plane stress. 
It will be shown that three types of problems may be distinguished. In certain special 
cases the stress distribution may be found independently of the velocity or thickness 
by solving three equations in as many unknowns. For general steady motion problems 
it will be necessary to solve six equations simultaneously for three stress components, 
two velocity components, and the thickness. These equations will be stated in Sec. 2, 
and reduced to a system of five first order, quasi-linear differential equations under the 
assumption of initial isotropy. For suitable boundary conditions, it will be possible to 
find the final stress, strain, and thickness distribution of the material directly, using 
the method of characteristics. The details of this method of solution will be described 
in Sees. 3 and 4. Finally, in See. 5, a step-by-step procedure for solving problems of 
non-steady motion will be briefly indicated. 

2. Basic equations. Let the sheet be referred to a set of Cartesian axes such that 
the x,y plane coincides with the middle surface, and let z = +}3h(2,y) be the equations 
of the bounding surfaces. Under the assumptions of generalized plane stress, the only 
non-vanishing averaged stress components are o, , ¢, , and 7,, . These components 
must satisfy the equations of equilibrium 


*Received March 29, 1950. 
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0 : 0 ‘ 
ay tee) a ay Te) = X, (la) 
0 0 ; 
art Ta) ay) = } (1b) 


Here X and Y are the components of body force in the zx and y directions, respectively. 
They may be functions of x, y, and h. 

In addition, if the material is to behave plastically, the stress components must 
satisfy a yield condition which may be written in the form 


rT) = 0, (1c) 


vy) zy 


Fi(c.,¢ 


Strictly speaking, Eq. 1(c) must be satisfied by the actual stress components. However, 
if the variation of the stress components across the thickness is not too great, it will be 
satisfied by the averaged stress components to a high degree of accuracy. 

In the particular case where fh is known, and where the boundary conditions in 
stresses are sufficient to set a fully piastic problem, Eqs. 1(a-c) may be solved directly 
for the stresses with no reference to the velocities." However, in general neither of these 
conditions will be satisfied, and it is necessary to consider the stress-strain relations. 
According to the “plastic potential” stress-strain law,” 


0u/Ox Ov/Oy 


Ou/dy + dv/dx 

2 a en, (1d) 
OF /dc, OF Oo, oF OF as 
Since F is a known function, Eqs. 1(d) provide two additional equations, but introduce 
two additional functions u and v, the averaged velocity components in the x and y 
directions respectively. The final necessary equation is obtained from the assumption 
that the plastic material is incompressible. For steady flow, this condition may be 
written in the form*® 


re) 0 
—(hu) + —(hv) = 0. (le) 
Ox OY 
Equations 1(a-e) are a set of six equations for the six unknown functions ¢, , ¢, , T., , 
h, u, and ». 
It will prove convenient to make the substitutions 
o, = 2k[w + x sin 26], Try = —2kx cos 26, 
(2a) 


a, = 2k[w — xsin 26], H =- log h. 
The function @ may be interpreted as the angle between a fixed direction and a principal 
direction (see, for instance, Ref. footnote 1), so that if the yield condition is isotropic, 
it may be written 

Foz, Cy, Tey) = x — J(w) = 0. (2b) 


1This has been done by the author in another paper. P. G. Hodge, Jr., Yield conditions in plane plastic 
stress, to be published in the J. Math. and Physics. 

2R. v. Mises, Mechanik der plastischen Formaenderung von Kristallen, ZS. angew. Math. Mech., 8, 
161-185 (1928). 

=R. Hill, Plastic distortion of non-uniform sheets, Phil. Mag. (7) 40, 971-983 (1949). 
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The substitution of Eqs. 2 into Eqs. 1 leads to the five equations 

(1 + f’ sin 26) w, —f’ cos 20w, + 2f(cos 26 6, + sin 26 6,) (3a) 
+ w+ fsin 20)H, — fcos20H, = X, 


— f’ cos 20, + (1 — f’ sin 20)w, + 2f(sin 26 6, — cos 26 6,) (3b) 
— fcos20H, + w — fsin 26)H, = Y, 

u,+tv, +uH, +H, = 0, (3¢) 

—2 cos 26u, + (f’ — sin 26)(u, + v,) = 0, (3d) 

(f’ + sin 26)(u, + v,) — 2 cos 26v, = 0, (3e) 


where subscripts now indicate differentiation. Equations 3(a-e) are a set of 5 quasi- 
linear first order equations for the 5 unknown functions w, 6, u, v, H. Once these quantities 
have been determined, the stress components are easily found from Eqs. 2. 

3. Characteristic equations. The characteristics of Eqs. 3 may be defined as those 
curves across Which it is possible for the derivatives of the unknown functions to exhibit 
finite jump discontinuities, the functions themselves being continuous. Let « = 2(s), 
y = y(s) be the parametric equations of a characteristic curve, and let w, 0, u, v and H 
be given along the curve as functions of the parameter s. Then, along the curve, 


u, dx + u, dy = du, v, dx + v, dy = dy, (4a) 


w, dx + w, dy = du, 6,dx + 6, dy = dé, 
(4b) 


H, dx + H, dy = dH. 
Equations 3 and 4 may be regarded as a set of 10 linear algebraic equations for the 10 
unknowns w, , w, , 6, , °-: , H,. In general they will possess a unique solution, so that 
a discontinuity in these unknown derivatives will not be permitted. An exception can 
occur only if the slope of the curve, dy/dx, is such that the determinant of the coefficients 
D, vanishes. By means of Laplace’s expansion,* the 10 X 10 determinant of D, may be 
reduced to the product of three determinants: 
Do = Dj: Ds’: Ds", 
where 


Di =|1+ f’sin20 —f’cos26 2f cos 20 2fsin 26 
—f’cos26 1-— f’sin26 2fsin 26 —2f cos 26 


dx dy 0 0 








0 0 dx dy ; 


dx dy|, 


‘See, for instance, A. C. Aitken, Determinants and matrices, Oliver and Boyd, London, 1939. 
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and 
pe” —2cos20 ff’ —sin26 ff’ — sin 260 0 | 
0 f’+sin20 f’+sin20 —2 cos 26 
dx dy 0 0 
0 0 dx dy 


Di, and Di” each vanish along either of the curves 
dy/dx = —cot a, dy/dx = tan £, (5a) 


where the substitutions 
i sin 2¥(|y| < w/yp), a=6+2, B=60-¢2 


have been made. Obviously these characteristics are real only if | f’| < 1; ‘n the re- 
mainder of this paper it will be assumed that | f’ | is actually less than one.” 
Each of the curves 5(a) represent a multiple characteristic. The final family of 


characteristics is obtained by setting Dj’ = 0: 
dy/dx = v/u. (5b) 


The curves 5(a) are called the first and second characteristics, respectively, and the 
curves 5(b) are called the streamlines. 

4. Conditions along the characteristics. The problem of finding conditions which 
must be satisfied along the characteristics is complicated by the existence of multiple 
characteristics. It will prove convenient to use one method for Eqs. 3(d, e) and find the 
velocity conditions along the characteristics. The velocity derivatives will then be 
regarded as known, and a different method applied to yield the stress relations along 
the characteristics and the relations along the streamlines. 

Consider the characteristic curves (Eqs. 5a) as curvilinear coordinates — = const., 


n = const., where £ and 7 are functions of x and y. It is easily shown that the relations 
n, Sina + n, cosa = 0, &. cos 6+ & sinB = 0, (6) 


must be valid along a first and second characteristic, respectively. 
The substitution of Eqs. 6, together with the relations u, = u;t, + u,n, , ete. into 


Eqs. 3(d, e) leads to the equations 


£. sin 8 Ue sina — Vv: COS a] ==. @ cos alu, cos 8 + v, sin 8] = 0, | 
(7) 
—&. cosa cos Blu; sina — v; cosa] — 7, sin Bsin alu, cos 8 + v, sin B] = 0. 


If Eqs. 7 are regarded as simultaneous linear equations for the expressions in square 
brackets, it is seen that the only solution is the trivial one. In view of the definitions of 


£ and 7, this may be written as 


du sin a — dv cosa = 0 (8a) 


‘The meaning of this restriction has been discussed by the author (see footnote 1) and by R. Hill 


(see footnote 3). 
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along a first characteristic, and 
du cos B + dv sin B = 0 (8b) 


along a second characteristic. 

To find the relations between w, 6, and H along the characteristics and streamlines, 
let us assume that wu and v and their derivatives are known from Eqs. 8(a, b) and con- 
sider Eqs. 3(a-c) as a set of three equations for the three unknowns w, 6, and H. The 
characteristics of these three equations will be those curves for which the determinant 
of the coefficients of w, , --- , H, vanishes. Let this determinant be denoted by A, , 
and let A, be the determinant formed by replacing the kth column of A, by the inhomo- 
geneous terms (i.e., the terms not containing w, , --- , H,) of Eqs. 3(a-c) and 4(b). 
In order for the derivatives w, , --- , H, to be able to exhibit finite discontinuities across 
the characteristics, it is necessary for all of the determinants A, , A, , --- , As to vanish. 
The first condition leads to the three curves given by Eqs. 5, of course, and the re- 
maining conditions lead to the relations 


(u cosa + vsina)[2f dé — cos 2p dw — (X sin B — Y cos 8)(dy/cos a)] 
= |[(u cos 8 + vsin B)w + (—usina + v cos a) f] dH (8c) 


+ (u, + v,)(f — w sin 2y)(dy/cos a) 
along a first characteristic, 


(usin 8 — v cos 8)[2f dé + cos 2b dw — (X cosa + Y sin a)(dy/sin 8)] 
= [(—usina + v cosa)w + (u cos 8B + vsin 8) f] dH (8d) 


+ (u, + v,)(f — w sin 2y)(dy/sin 8) 
along a second characteristic, and 
(u. + v,) dy+uvdH =0 (Se) 


along a streamline. 

5. Conclusion. Equations 5 and 8 may be replaced by finite difference equations and 
used to obtain numerical solutions of boundary value problems. There are a great num- 
ber of possible such problems, corresponding to various combinations of the Cauchy 
and Riemann problems for a system of two hyperbolic equations. The details of the 
method may be readily adopted from corresponding treatments in compressible flow.° 

For the general case of non-steady motion, it is necessary to first solve the problem 
at the instant ¢) when plastic flow begins. Since H is an initially known function, this 
involves solving Eqs. 3(a, b, d, e) for w, 6, u, and v. This problem is similar to the one 
discussed in the preceding sections, and, in fact, (5a) and 8(a, b, ce, d) may be used, 
since they are not dependent upon the motion being steady. The right hand side of 

‘See, for instance, W. Tollmien, Stationaere ebene und rotationssymmetrische Uberschallstroemungen, 
Technische Hochschule Dresden, 1940. This work has been translated by the Air Materiel Command, 
Wright Field, Dayton, Ohio as translation A9-T-1, Contract W33-038-ac-15004 (16351), Brown Uni- 


versity, Providence, Rhode Island. 
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Eq. 3(c) must now be replaced by —/H, , and, since the left side of Eq. 3(c) is known at 
t = t) , H,(t)) is known throughout the field. A finite difference approximation then 
yields H(t, + At) throughout the field, so that the process may be repeated. Thus a 
step-by-step solution may be obtained for all values of ¢. 


A NOTE ON THE PRINCIPAL FREQUENCY OF A TRIANGULAR MEMBRANE* 
By G. POLYA (Stanford University) 


The principal frequency of a membrane of triangular shape is exactly known in 
two simple cases: for the 45°, 45°, 90° and the 60°, 60°, 60° triangles.’ As will be shown 
in this note, an exact solution of comparable simplicity exists also for the 30°, 60°, 90° 
triangle—a result which, to the author’s knowledge, has not been observed before. 

The three lines the equations of which in rectangular coordinates x, y are 


y = 0, z = a3'’’/2, y = 13°”? (1) 


delimit a triangle, one half of an equilateral triangle with side a. Define 


u = sin ry — 37'2x)/a - sin x(y + 37'2)/a - sin 2437'?2/a 
(2) 
- sin r(y — 3'”x)/3a - sin r(y + 3'”’x)/3a - sin 2ry/3a 
By elementary transformations we find that 
—32u = cos 2n(5y + 3'x)/3a — cos 2r(4y + 3'?2x)/3a + cos 2x(y + 3'°3x)/3a 
(3) 


— cos 2r(5y — 3'x)/3a + cos 2n(4y — 3'’2x)/3a — cos 2x(y — 3'’°3x)/3a. 
We see from (3) that u satisfies the equation 
Ure + Uy + 9 '112r'a*u = 0, 


from (2) or (3) that uw vanishes along the lines (1), and from (2) that u does not vanish 
in the interior of the triangle between the lines (1). Therefore, u represents the principal 
mode. 

The well known solution for the two other triangles mentioned at the beginning 
can be presented in a strictly analogous form. Any of these three triangles, repeated by 
successive reflections, covers the whole plane without overlapping and there are no 
further triangles of this kind. Therefore, there are no other triangles for which the 
solution can be presented in a comparably simple form as a function of x, y without 
singularity in the whole plane. The last remarks indicate also the heuristic reasoning 
which led to the solution (2). 

*Received April 17, 1950. 
1Sponsored by the Office of Naval Research. 
*Rayleigh, The theory of sound, Dover Publications, New York, 1945, §199. 
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THE RANKINE GAS FLOW IN THE HODOGRAPH PLANE’ 
By MAX M. MUNK (White Oak, Maryland) 


Introduction. In the theory of plane two-dimensional potential flows of a perfect 
incompressible fluid, the flow associated with the superposition of one point source 
and of the parallel flow is known as the Rankine flow. When transforming into hodograph 
variables and when choosing the Legendre reciprocal potential L as the dependent 
variable, the Rankine flow is represented by a sole point source positioned at a velocity 
point representing the velocity at infinity in the physical plane. 

It is the aim of this paper to set up a corresponding hodograph solution for a perfect 
gas, arriving with subsonic velocity and flowing past the head of an obstacle extending 
infinitely far downstream. We shall have to employ the stream function y in the physical 
plane which is associated with the physical potential ¢, and we shall have to employ 
and will primarily deal with the Legendre reciprocal potential L in the hodograph plane, 
which is associated with a reciprocal stream function &. The functions ¢ and é will 
not be used in the numerical computations. | 

We recall that the sets of differential equations for these quantities in terms of the 
hodograph variables w and @ are the following ones, with Q being written for 


Q= (1 - a w/a — w) 


by way of abbreviation: 


10 _ ak 
wae” aw’ 
(1) 
0 _ _Qp al 
dwt sO” 
1 dd _ 10y 
wae p dw’ 
(2) 
a _ _Qdy 
ow p 00° 


The velocity w is expressed in units consistent with the following. The velocity incli- 
nation angle is denoted by @ (its square is often denoted by 7), p denotes the density 
p = (1 — w’)’’"~” and y denotes the expansion exponent of the gas. The Legendre 
transformation employed relates to the Cartesian coordinates in the physical plane and 
in the hodograph. It does not relate to polar coordinates. 

Equations (1) and (2) possess particular solutions of the form 


L(w, 0) = L,(w) exp in8, 


E(w, 0) &,(w) exp i(n@ + 2/2), 


Y(w, 0) = y,(w) exp i(n@ — 2/2), 


*Received Apr. 19, 1950. This work was sponsored by the Office of Naval Research. 
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L,, being a solution of 


w s | ow dl — n pQL, = 0 (la) 


dw du 


In this paper, only integer values of the parameter 7 will occur. The special solution 
of (la) which becomes zero at w = O will be designated as L, . Replacing in such a 
solution the integer n by —n does not lead to a new independent solution. The symbol 
L_, will be used in this paper to denote a second solution of (la) which is independent 
of the first solution LZ, but which is not its analytic continuation with reference to n, 
and which is not obtained by substituting —n for +-n. In this paper, n is always positive. 
Substituting for the subsonic velocity w a “distorted” velocity s defined by 


ds/s = ~/Q dw/w, (4) 


We transform the two sets (1) and (2) into similar sets (1b) and (2b) superior to (la) 
and (2a) in symmetry, in that then in each pair of equations the same factor appears 


on the right hand side: 


1 dé — dL 
sag ~ °V 2 a5 
‘3 (16) 
OF _ eVQ aL 
dss 0@’ 
l 0g l 4 0yY 
sn Vie’ 
(2b) 
0d ms /Q oy 
ds 98 00° 


Equations (1b) and (2b) can be interpreted as relating to the steady flow say of 
electricity or of heat in a disc occupying the distorted hodograph plane. The disc has 
‘tivity (or resistance) the conductivity being a function of the distance 
from the center but not of the flow direction. This conductivity is proportional to 
pv/Q for (1b) or to ~/Q/p for (2b), and becomes zero at the sonic speed Q = 0. 

This indicates that for small variations of s, the solutions in the distorted subsonic 
portion of the hodograph must approach solutions of Laplace’s equation. If the re- 
sistance changes are very gradual compared with the changes of the potential, an ap- 
lution can be obtained from a solution of Laplace’s equation by increasing 
the potential gradients and decreasing the stream densities in the same ratio given by 

f the conductivity. This thought leads to the so-called asymptotic 


the square root 
particular solutions of (1) and (2) which are approached by the exact solutions as n 


approaches infinity. The asymptotic solutions are 


variable condi 


proximate s 


L = p Q'S” exp in6 (5) 
and similarly for the stream function. In the supersonic regime an analogous relation 
holds, the wave equation replacing Laplace’s equation, and the sine or cosine of the 
logarithm of s replacing the power of s. These asymptotic relations break down at the 


sonic speed. They can be obtained or verified by inserting the multiplier of the exponential 
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function in (5) into (la), and setting the terms containing the two highest powers of n 
equal to zero. This gives two ordinary differential equations for s" and its factor. 

The asymptotic behaviour of the reciprocal potential L and of its associated stream 
function y clarifies the character of the infinite series in terms of L, or ¥, to be established. 
Their convergence and the type of singularity represented is determined by the terms 
of higher order and just these terms are described by the asymptotic expressions. The 
series are thus recognized to behave in many respects like the corresponding series 
relating to Laplace’s equation. 

Sets (1) and (2) are closely related to each other and in consequence the functions 
L, and y, are also closely related. The following relations hold for all particular solu- 
tions, whether designated by n or by —n (n > 1): 


w OL, 
vn = o| nb, ~~ n Ow | (6) 
. 
os w av, 
n= p(n” — 1) ny, , n Ow |’ (7) 
y en E on 
SO) .. pan ~ 1 Oe. (8) 
Ow On 


Note that for the first solution the stream function multiplier associated with the 
reciprocal potential is 
L. a 2 
i. = Bi _ — 1)w"F(a, b;c; w’) (9) 
where the hypergeometric series is set up with parameters a, b, and c which are approp- 
riate for the stream function and for the value of n. 

In (6) to (8) n has to be taken positive if so taken in (3). 

The Rankine flow. In the theory of the potential flows of the ideal incompressible 
fluid we are attracted by the Rankine flow on account of its simplicity from the mathe- 
matical view point. In the complex plane z = x + 7 y representing the physical plane 
the potential is the real part of a complex variable P 


P = logz+2z 
¢ = dP/dz (10) 


L —-P+2¢ = log(l- O+1 


Equations (10) presume unit velocity at infinity and unit source velocity at unit distance 
from the source. The source is at the origin of the coordinate system x, y. The impact 


point is at x = —1, the flow is in the x direction. The vertical x = —} is the line of 
constant velocity magnitude equal to the velocity magnitude at infinity. In the hodo- 
graph this line is transformed into the circle w = 1, containing the hodograph source. 


This circle will be called the singularity circle. 
We have to elect the features of the flow of the incompressible fluid which we wish 
it to share with the corresponding gas flow. We are guided in this choice by expediency 
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with regard to the mathematical treatment. We emphasize the condition at infinity. 
We not only require there the velocity to be constant and parallel. We furthermore 
require that the conditions at infinity be represented in the hodograph plane by a single 
source. We find it finally desirable to derive the series for Z in terms of L, directly from 
the corresponding series for the incompressible fluid. It is therefore desirable to require 
that at the singularity circle a suitable series should agree term by term with the corre- 
sponding series for the incompressible fluid. A quantity acceptable for that is dL/dw. 
This is equivalent to the requirement that the lines connecting pairs of points of equal 
velocity-direction and of undisturbed velocity magnitude relating to the incompressible 
fluid and to the gas be normal to the velocity in question. This implies that the curve 
of undisturbed velocity magnitude for the gas have a vertical asymptote x = const. 

The solution for L is given by two series applicable respectively to the inside and 
the outside of the singularity circle. They are determined (a) by the condition that each 
term be a solution of (1), (b) that the differential quotient with respect to w of each 
term be equal to the corresponding quantity applying to the incompressible flow, and 
(c) that the terms of the two series approach the same value at the singularity circle 
except at the singularity point itself. Let the ‘‘prime” denote differentiation with respect 
to w. Then for the incompressible flow we have 








5, = — >» (1/n)w" cos n6 
n=1 
(11) 
Si = — >> w"" cosné 
n=1 
S, = log w — >> (1/n)w™ cos n0 
n=1 
(12) 
Si; = >> w"" cos né 
The corresponding series for the gas flow are 
L,, cos n6 
8, = —2,- 7? 
, > [Lal 
[L0/,+ LL, 1 | 
=- in To. — | L, cos nd 
a > Fava ~~ ee 
(14) 


7 2h, | 
+ pH Fee i 47 |b cos n@é 


The subscript s denotes the value of the quantity at the singularity circle, which 
latter is subsonic. 

Concluding remarks. The position in the physical plane associated with any hodo- 
graph point is equal to the gradient of L. The stream function y is computed by applying 
(6) and (9) to the solutions (13) and (14) term by term. 

Modifications of the solution can be obtained by adding through the entire hodograph 
plane a finite number of terms containing the L, (n positive) or even an infinite series 
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of such terms converging appropriately. The method can also be applied to other singu- 
larities of integer order. 


APPENDIX. On the choice of the second solution if n is integer. 

The differential equation for L,, is related to the hypergeometric equation. It imposes 
a definite ratio on pairs of consecutive terms of a power series representing a solution. 
This ratio depends on the coefficients of the differential equation and on the order or 
power of the two terms. It is therefore a straight-forward procedure to write down a 
series P,, extending infinitely to the ascending and to the descending side complying 
with this requirement, starting with a term Ax’ where A and @ are arbitrarily chosen. 
Any power series actually representing a solution must be a special case of R, , namely 
o must have the proper value. This special value is computed from the so-called indicial 
equation, obtained by setting the multiplier of any one term of R, , preferably the term 
of power ¢ — 1 equal to zero. The series breaks then off on one side, in this case on the 
descending side, and the terms in ascending powers may form a convergent series. The 
indicial equation is of second degree and ordinarily furnishes two roots, namely ¢ = n 
and ¢ = —n in the present case. 

Differentiating R, with respect to o leads to a second series (not a power series) 
which under proper circumstances may also represent a solution. 

If n is an integer, the two series derived from R,, have terms of equal power in com- 
mon, and therefore at best would represent substantially the same solution. In this 
case, the series 0R,/d0 furnishes a second solution independent of the first solution. 
But any product R, multiplied by any function of ¢ would likewise, and when differ- 
entiated would not furnish the same second solution. The second solution is general, in 
that it contains two integration constants. 

No substantial reasons have become known to this writer why to select one particular 
second solution in preference to another. The solution ordinarily considered as standard 
is derived from R,, in the manner described and not from R,, multiplied by a function of 
c. In the series representation of that solution, the hypergeometric series appears 


Pp 


multiplied by log (A 2”), where A = 1. 

NoTE ON THE CONVERGENCE OF (13) AND (14). By A. Van Tuy. These series are 
convergent as a direct consequence of the asymptotic expressions for L., and L1, im- 
plied in (5). The proof of these asymptotic relations follows from a more general result 
of J. Horn [4] from which we have 


| - reget + o(2)|, 
n 
Lo. = ns gel 1 + o(2) |. 
n 
Denoting the value of s at the singularity circle by s) , we find 


fe (et). [eee Ade = Che) 
i Nes a Le. Lv whe \) Ae 
org |en= 
F > | = NE) Nal 
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From this follows the convergence of (13) for s < s) , the interior of the singularity circle, 
and the convergence of (14) for s > s) . The reader may wish to compare this with 
reference [1], pp. 46-49. 
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Ordnung mit einem willkuerlichen Parameter, 


ON THE METHOD OF INVERSION IN THE TWO-DIMENSIONAL 
THEORY OF ELASTICITY* 


By E. STERNBERG anp R. A. EUBANKS (dllinois Institute of Technoloqy 


1. Introduction. The method of inversion, originally introduced by J. H. Michell [1], 
has led to a variety of technically significant solutions to “plane” problems in the 


theory of elasticity [2], [3], [4], [5]. The usefulness of Michell’s stress-field transformation 


stems from its invariant properties which assure the preservation of an important class 


of boundary conditions. In the present note we show that any conformal stress-field 


transformation which preserves the principal-stress trajectories for every choice of the 
antecedent Airy function, is essentially a Michell transformation. 

2. The Michell transformation. The inversion theorem of Michell may be stated as 
eal and biharmonic in a rezion R of the z-plane, and let R* be 


follows. Let U'(z, Z) be 1 
the image of R with respect to the mapping’ 
az + b a b 
¢ = w(z) = <a = |. (1) 
cer d Cc d 
Then the function 
U*(g, F) = AU[G(O), gO), (2 


where hh” = | w’ |” and gq is the inverse of w, is biharmonic in R*. The stress fields gene- 


rated by U and U*, considered as Airy functions in R and #&* respectively, are related 
according to ~ 


o* + ir* = No + ir) +2, 


l ‘ ‘ 
\= p = ru — MU, — 2.02). | 


*Received Mar. 15, 1950. 


\Michell actually used ¢ = 1/2. 
*Subscripts attached to functions which originally bear no subscripts denote partial differentiation 
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Here [c¢, 7] and [o*, 7*] are the normal and shearing stresses on any are I and on its 
image I'* with respect to the mapping [ = w(z). 
Moreover, the stress-field transformation characterized by the mapping (1) and the 
law of Airy function association (2) has the invariant properties: 
(A) The images of the principal stress trajectories of R are the principal stress 
trajectories of R*. 
(B) If any are T of R is acted on by constant normal tractions only, so is its image 
oe. 
(C) A concentrated load acting at a point z, of a boundary arc I of R, and including 
a certain angle with I, is carried into a concentrated load of the same magnitude 
acting at w(z,) and including the same angle with the image are I* of R*. 
Properties (A) and (B) were also established by V. P. Jensen and D. L. Holl [6] 
by aid of derivatives of non-analytie (polygenic) functions [7], [8]. 
3. A converse of the inversion theorem. We now prove the following theorem. Let 
R be a region of the z-plane and let R* be the image of R with respect to the conformal 
mapping 
¢ = w(2), w'(z) = “ = he’ #0. (4) 
Moreover, for every Airy function U(z, z), bi-harmonic in R, let there exist an Airy 
function U*(¢, ¢), bi-harmonic in R*, such that the corresponding stress-field trans- 
formation preserves principal stress trajectories, and the image field of stress is purely 
hydrostatic only if the antecedent field has the same property.* Then w(z) is given by 
(1) and U* is given by (2), i.e., the transformation is a Michell transformation. 
To establish the theorem, we recall a result of Jensen and Holl [6] who showed that 


o + ir = yulz, 2, 6) 
(5) 
= 2U..+ wun", ) 
where 
Hos.s) = 2U; (6) 


and y, is the directional derivative of H along I, 6 being the inclination of I’. In view 
of (5), property (A) is equivalent to the statement 


vu+[w(z), w(z), 0 + 6] - Yn (7) 
whenever _ 
Yu(2, a 6) = Yu ’ (8) 
provided - 
H*(t, ¢) = 2U?. (9) 


Equations (7), (8) by aid of (5) become* 


y ry 440 
U un = l 22€ ’ 
(10) 
“a Tx 44(6+5) 
Uz; = Ulex . 
3This restriction is essential in order to rule out the trivial transformation which carries all antecedent 


stress distributions into hydrostatic fields of stress, 
‘The partial derivatives of U* with respect to ¢ and ¢ are to be evaluated at ¢ = w(z) 
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Thus 
U,,U?(w’) = UnUt(w’)’, (11) 
which implies that the function 
o(z, 2) = Us, U}(w’)’ (12) 
is real-valued. Furthermore, assuming that the original stress distribution is not hydro- 


static, so that U., * 0, it follows by hypothesis that U}, ~ 0 and hence @ does not 
vanish identically. For convenience, let 


gz, 2) = UU a/. (13) 

Equation (12) then appears as 
U,, = ow’) U3, , (14) 
where ¢ is again real-valued. Equation (14) constitutes a necessary and sufficient condi- 


tion for the preservation of the principal-stress trajectories. 
We next apply to (14) the condition that U(z, z) and U*(¢, £) are both biharmonic, 


i.e., U,.;; = 0 and U*,;; = 0. This leads to 
w'! p w'! (Qv.)° 7 - 
P i + 20. | l 22 + E — 2 cea Bel] l oo @= 0, (15) 
L” w w ro) 


or, since (15) must hold for every bi-harmonic U, 


ge 7 + 2%, = 0, (16) 
w'! o . o-~ 
Wa OT - Me = 0. (17) 


The complete solution of (16), subject to the requirement ¢ = @, is 


(z, 2) = x(w’w’)'” = K/h (18) 


A 
Ye \e; 


with « an arbitrary real constant. Noting that (16), (17) require ¢,, = 0, we conclude 


from (18) that 
wl! , 1 ta] 
. — oie _ QO) 
(2 2 \w’ ) 0, ai 
i.e., the Schwarzian derivative of w(z) vanishes. The complete solution of (19) is given by 


az +b 

w(z) = = x 

@) cz + d’ (20) 

and since the mapping is to be (1, 1) we may put ad — be = 1. In order to arrive at the 

law of Airy-function association we now integrate (14) by use of (18). This integration 
yields, 


ue =* + Ast + od + Ge + Bl, (21) 
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where A, B are arbitrary real constants, and a@ is an arbitrary complex number. It is 
readily confirmed by direct computation that the constants A, B, and a give rise to 
an arbitrary, uniform hydrostatic stress-field in the image domain R*. The constant x, 
on the other hand, affects merely the scale of the image stress-distribution. We may 
therefore put A = B = a = 0, x = 1. This completes the proof. 
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A MINIMUM PRINCIPLE FOR STRUCTURAL STABILITY* 
By H. J. WEISS ann G. H. HANDELMAN (Carnegie Institute of Technology) 


1. Statement of the problem. In a recent paper,, W. Prager has discussed the 
problem of structural stability from the following point of view. Consider a given con- 
figuration of a deformable body, referred to a fixed system of rectangular axes, 
z,(i = 1, 2, 3), under a set of stresses \o,; which are in equilibrium with given surface 
tractions. These stresses are prescribed only to within the arbitrary constant factor X. 
The configuration is assumed to be stable if \ is sufficiently small. A system of in- 
finitesimal perturbation displacements u, is then applied, and the question is asked for 
what values of the factor \ will the equilibrium become indifferent. 

The solution to the problem leads to the following system of linear, homogeneous, 
second order, partial differential equations 


Jiic = 0, (1) 
subject to the homogeneous boundary conditions on the surface 

Jin, = 0, (2) 
where 


J ii = [r<3 + BA (Cine; — F ip€ps) = AG py i]. 





*Received May 29, 1950. 
1W. Prager, The general variational principle of the theory of structural stability, Q. Appl. Math., 4, 
378-384 (1947). 
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The quantities appearing in the definition of J,;; are defined below: the strain, «;; , 
associated with the perturbations wu, is given by 


é:;; = 4(u,,; + U;.;)3 (3a) 


the additional stresses produced by these strains satisfy the linear relation 

3 eo (3b) 
in which C;,,;, 
the local rotation due to the displacement u, , 


is symmetric with respect to 7 and j and with respect to k and 1; w,; is 


Wis = $(Us,5 — U;,0)3 (3c) 
and n,; is the unit vector normal to the original equilibrium surface. The problem dis- 
cussed here will be restricted to the elastic case; thus 


\ 


y ory a ¢ V ° ~ 
Cin = 26,( 6.45; = — Sichusd, 


ad 1 — 2p 
where G, is the modulus of rigidity and v is Poisson’s ratio. It should be noted that 
Cx. is symmetric with respect to (77) and (kl). The initial stresses o;; are symmetric 
with respect to 7 and j and satisfy the equilibrium equations 


O:;,4 = Y. 


Prager has shown that the differential equations (1), subject to boundary conditions 
(2), are equivalent to the statement that the first variation of the integral 


| [C sipe€ii€pa TH AGis(Up, Up.; — Ei€pi)] dV (4) 
| 
ranishes. Integration is carried out over the volume V of the original unperturbed body 
and only the displacements u, are varied. 

2. The minimum principle. It will be shown that the eigenvalues \ and the corre- 
sponding eigenfunctions of Eqs. (1) may be obtained from a minimum principle applied 
to an integral closely associated with that of Eq. (4). Since the problem in question is 
one of stability, only the first eigenvalue will be considered. 

It is more convenient to express Eqs. (1) and (2) in terms of the displacement u, 
directly. Thus, from Eqs. (3a), (8b), and (3c), J;; may be written as 


r 
/ ‘ P ‘“ a 
[o,:(3u;,, — Uy, 3) — Op;(Ui.p + ae (5) 


J ij =C pallp.g + 4 L 


and the integral (4) becomes 


p " : 
[ ‘iva, Uy.a T 4 a (Buy, Mp, j — Uji ,gUj,p — Up, Uj,p — Ui,pUly, i) dV. (6) 


Diu, r) = | Cyipdts.Yy.2 AV, 


Vv 
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1 f ' 
H(u) -, 4 7; (BUy, Up, § — Ui, Qj ,p ~~ Up, itp Us ,pUy, i) dV, (8a) 
e V 
\ l 7 
H(u, v) = “s 6;5;(3Uy, Vp. — Us.Yinp — Up. Vico — Ui.) AV, (Sb) 
4Jy 


it is seen that 


Diu + v) = | Cy inating + U:,;)(Up.e H Yp.q) AV 


Jv 


Diu) + Div) + | Cying(Ui.Up.g H Up. ats.;) AV 


/¢y 


II 


= Diu) + Dv) + | IC} spetti. Me.e H Coatits.®%r.e] AV, 


\ 


Diu + v) = Diu) + 2D(u, v) + Div); (9) 
and 


. 


l , 
H(u +. H(u) + Hv) - | 6; |3(Uy, Up. HU p.iUp.3) — (Ui pip HF Vie, p) 


— (Uy. Vj .p Hepp H Ui plo.s HUspty.;)] dV, 


Hu +i H(u) + 2H(u, v) + Hv), (10) 


since o,; is symmetric. 

The vector function u; = u,(2, , % , 2) Will be called an admissible vector if it is 
continuous and has continuous partial derivatives of at least the second order. The 
first eigenvalue of the differential equation (1) subject to boundary conditions (2) can 
then be found” from the following minimum principle: 

Of all admissible vectors u,; , that one which minimizes D(w) under the side 
condition H(u) = 1 is an eigenfunction of the differential equation (1) 
subject to boundary conditions (2). This minimum value of D(u) is the 
corresponding eigenvalue X. 
It should be noted that this minimum principle can be rep aced by one requiring that 
D(u)/H(u) be a minimum. The corresponding eigenfunction will then be determined 
only to within a factor of proportionality. 

In order to prove this statement, let u, be the minimizing function and X the corre- 

sponding minimum value. Then 


D(u) = AH (u). (11) 


Let £; be an arbitrary admissible vector and ¢ an arbitrary parameter. Then, since u, 
minimizes D(u)/H(u), it follows that 


Diu + e&) = AH(u + é€). 


2R. Courant and D. Hilbert, Methoden der mathematischen Physik, vol. 1, Julius Springer, Berlin, 
1931, pp. 345-348 
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From Eqs. (9), (10), (11), this can be rewritten as 
as hs ek 
2e) Diu, §) — AH(u, §) + 5 [D&) — AMN()] 7 2 O. 


Since ¢ is chosen arbitrarily, the above relation will be valid only if 


D(u, &) — »H(u, &) = 0. sie 
Now 
Du, &) — AH (u, €) 
- i lc. palli.ibo.g 1b | F6j(BUp,sbo,5 — Usixbi.n — Un.ikinn tabs.) | oe 
<2 fe ee oe . 7 
” sii ty 0; ;(BUy, &.5 — Ui.rti.p — Up.itinp — i) ae ae 


By Green’s theorem, this can be transformed into 


pour 
r . 
1C; pall, af nj; + 4 [o;;(3u, i Uy pen; sid o 5 ;(U; Dp a i». dma] ds 


7s \ 


ne | 1C; palyp,a E; + : [o,;(3u,,< — U; »)]. i€> _ [o, (Ui.p + vs d]L8) dV, 
Jy \ 


where S is the bounding surface of V. By an appropriate interchange of subscripts, 
Eq. (12) can be written finally as 


‘ X 
| | Cana q + * (Bs. 7 Uy, i) —_ 4 Oyi(Up, j + Uu; pus dS 


/Ss 


a | {c pally,g , [o,;(3u;,, = Oe Tyi(Uy,; + ws ab dV = 0. 


jy (ern TE 
Now the vector &, is chosen arbitrarily in the volume V and on the boundary S; conse- 
quently 

f r 

\Cisvatt ‘ 7. 4 [o,;(3u;,, aid Uy, i) ery Tyi(Uy, j + U;.»)] ss 0, 
and 


. 
‘Coin q a 4 [o,;(3u;,» = U, i) rs Op i(Uy, ; + Wp = 0.. 


These are identical with the differential equations (1) and boundary conditions (2). 
Thus the solution of the minimum problem yields an eigenfunction and corre- 
sponding eigenvalue of the differential equation. It can also be seen that if an admissible 
function v, is used to generate an approximate solution, D(v)/H(v) =  , the actual 
eigenvalue. 
3. Quasi-orthogonality and reality of the eigenfunctions. A quasi-orthogonality 
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relationship between the eigenfunctions of Eq. (1), from which the reality of the eigen- 
values follows, can be established by means of a standard technique.* Define 


M j;(u) = Ciiply.e 
and 


N,,(u) _ —3(30,,u;,» — Gigly,5 — FjgQki,p — O jglly, i) 


Then 
J 53 = M,,(u) asl AN, ;(u). (13) 


Equation (1) takes the form 

[M;,(u) — AN.,@)],, = 0, (14a) 
and boundary conditions (2) become 

[M,;(u) — AN,;(u)]n, = 0. (14b) 


Let u,; and v, be two eigenfunctions and \ and yu their respective eigenvalues, where 
\ * uw. Then 


[M.,(u) — AN,,()]., = 0 


and 
[M.,(v) — uN,,(v)],, = 0. 


Multiplying the first by »; and the second by u; , subtracting, and integrating yield 
[ Jude) — uyJSa.dde¥ = 0. (15) 
This can be rewritten by Green’s theorem as 
J toidala) = wi ses dS = fey. Fale) — wp] dV = 0. 
But J,;,;n; = 0 on S; therefore, the surface integral vanishes; also 


| [v;, .M,;(u) “a u;,.M,;(v)] dV _ I [v;, Ciivatty. oo U; iC sinWr.al dV 
Vv 


V 
aad I Csivald;. Up. _ Uy. a¥i,s] dV =0. 
Vv 
Therefore, Eq. (15) becomes 


r / v;, «NV, ;(u) dV = uf u; iN. ;(v) dV = 0. 
V Vv 





3L. Collatz, Eigenwertprobleme, Chelsea Publishing Co., New York, 1948, p. 61-65. 
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Now 
, ; ¢ ' 
v,:N,,(u)dV = — I v (36 ;,U;,.p — Ciglp,j — Cjpi,p — Fipty,:) AV, 
V vy 
1 f , * 
= — | U;, (SOV; .p — Figy.5 — Cighi.s T;Wy.i) AV, 
Jy 
= | u;,.N,,(v) dV. 
V 
Therefore, 


(A —p) | 0,,.Niw dV = 0. 


Since \ =~ u, it follows that 


4 


I v,;,.Ns(u) dV = 0. 
Vv 
Now let us consider 


| v;[M,,(u) — AN, ;(u)],; dV =0, 


b 


which can be recast as 


. . 


| v,[M,;(u) — AN,,(W]n; dS — | v;,,M;,,(u dV +A | vy, .NV;;,(u) dV = 0. 
' Jy 


Ss 


But 
[M,(u) — AN;,(w)]n; = 0 


on S and 


D oaNulu) AV = 0. ™ 
‘ 
Therefore, 
[ vty av = 0. - 
vy 


Equations (16) and (17) are the quasi-orthogonality relations which hold between the 
eigenfunctions. It should be noted that these results are somewhat more complicated 
than the weighted orthogonality conditions usually found. This depends on the fact 
that the eigenvalues appear in the original differential equations and boundary condi- 


tions as coefficients of terms containing the highest order derivatives of the unknown 


functions. 
Now assume that A is complex; that is, \ = s + 7t where s and ¢ are real, and ¢t ¥ 0. 
Corresponding to \ there will be a complex eigenfunction u; = o; + i7; where o,; and 


7; are real vector functions and 7; # 0. Since Eqs. (1) and (2) are linear with real co- 
efficients, it follows’that \ = s — it is also an eigenvalue with corresponding eigen- 
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function u; = o; — ir; . A computation, similar to that used in the development of the 
quasi-orthogonality relation, will show that 


(X—2) | wu dV = 0. 


Since N,;(u) is linear in the derivatives of u; , 
N,,@) = Nie — ir) = N,,(0) — iN;,(7) 


and 


[ (o;,, + i7;,)Ni(o — ir) dV = [ [o; Ni;(o) + 7;,.:Ni;(7)] dV 


V 


~i ff foNult) — raul] av. 
The imaginary part of this integral vanishes by the quasi-orthogonality relationship, 
and the integral reduces to 


“ 


I o,..Nilo) + 1y..Nild] av. (18) 
“Vv 
If this integral does not vanish, then \ — \ = —2it = O and ) is real. Thus the eigen- 
values and eigenfunctions will be real if the integral appearing in (18) is non-zero. 
Now 


| T,:N;;(7) dV = an | 0:;(37.;Tp.i — Ti.pTr.i — Tp.iTi.yp — Ti.vTi.2) 2V = H(7), 


and 
| o;..Ni(o) dV = Ho). 
dv 


According to the minimum principle established in Section 2, 


\ = min |D(u)/H(u)). 


It has been assumed, tacitly, that a solution exists to this minimum problem. Under 
these conditions, H(u) = 0 only if D(w) = 0. But D(u) is the strain energy of the per- 
turbed state and will vanish only when the perturbation is a rigid body displacement. 
Consequently, H(u) # 0 if the minimum problem is to have a finite solution. Similarly, 
since o, and 7, are eigenfunctions of the differential equation (1), neither H(o¢) nor 
H(r) can vanish if there is to be a finite, non-trivial solution. Integral (18) will then 
vanish only if /(¢) and H(r) are of opposite signs. If H(u) is of one sign for all eigen- 
functions, u, , the eigenvalues and eigenfunctions must therefore be real. This restriction 
is not artificial, however, for H(u) is the difference’ between the work of the surface 
stresses and the work of the body stresses \o,; on the strains associated with the displace- 
ments u,; . Thus if this work is always positive or always negative, the eigenvalues and 
eigenfunctions will be real. It is also seen that if H(u) is of one sign for all admissible 
functions, the eigenvalue \ will be of the same sign since D(u) 2 0. 


‘This relation was pointed out to the authors by Professor W. Prager. 
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NOTE ON AERODYNAMIC HEATING* 


By H. W. EMMONS (Harvard University) 


The acceleration of an initially stationary semifinite viscous fluid by the sudden 
starting from rest of the infinite bounding plane at velocity U parallel to the plane has 
been derived many times. See for example, Aerodynamic Theory.’ The equations of 
motion for this case assuming constant physical properties reduce to the well known heat 
conduction equation, 


Ou ou 
at = | ay” ’ (1) 


yielding as the solution for the velocity distribution 


uz U 3. [ exp [—6”] dB. (2) 


us “y/2(vt)*/* 


The question of the temperature distribution resulting from the energy dissipated 
has not been previously discussed to the author’s knowledge and seems worthy of a note. 
The energy equation reduces for the present case to 
or _ a 7.2 (2) (3) 
at dy’ ° c, \dy/’ 
where the physical constants a, v, the thermal diffusivity and kinematic viscosity, show 
the rate of spread of the energy and the momentum respectively. Their ratio is the Prandtl 
Number 


Pr= =". (4) 


A particular integral 7, of (3) may be written down by the method of sources 




















twol T tt « — "I ___drdrh__ 
1 i I. Cy (24) exp { 4a(t — r)) [4ra(t — 7)]'” (5) 
On substituting the velocity distribution from (2) 
ete it Tf fy = »)" _ ao Sw 
= TCy i T l. exp { 4a(t — 7) Qvr} [4ra(t — 7)]'” (6) 
This result may be simplified as follows: 
Let £ = —(y — d)/[4a(t — 7)]'” be used to replace A. One integration then yields 
U? rt { y? v Dy I’ ” 
1, = TC, | xP \~ 4at — 4az + 2Qyr)_7r(4at — 4ar + 2y7) dr. (7) 
Now let ¢ = y’/(4at — 4ar + 2vr7) 
vy? (2/Pr)(y*/4at) e dt 
T, = pr/(2 - Pry”? [ , 8 
xc, | a ¢[(4atz/y*) — 1)" ®) 


*Received June 2, 1950. 
IW. F. Durand, Aerodynamic Theory, vol. III, p. 64. 
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/ 


A trigonometric substitution ¢ = y’/4at sin’ 6 gives finally 


7? ‘ /2 : 
T, = — [Pr/(2 — Pr]” : exp {—y’/4at sin’ 6} dé. (9) 
Cp W Jgin-*(Pr/2}*/* 
Since there were no sources except those of dissipation, this becomes the solution for the 
adiabatic plate if we add Eq. (9) to the initial temperature T, . We take then as the 
particular integral the solution for an adiabatic plate 


72 w/2 
T,=7,+ 2 [pr/2 — Pry” 2/ exp {—y?/4at sin? 6} d@. (10) 
Cp TW Jein-*(Pr/2}*/* 


To satisfy any other boundary condition at the plate, it is only necessary to add the 
appropriate complementary function chosen from the many known solutions of the heat 
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conduction equation. In particular, for a rate of heat transfer g we add the complementary 
function 


7. = 22 (at/x}'? exp {-y'/4at} — Ye [exp [-"] a8). (A) 


w/2lat}*/* 


Returning now to a consideration of the solution for an adiabatic plate, we note that 
while the integration must in general be carried out numerically, the most important 
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question, that of the plate temperature, can be answered without difficulty. If y = 0 
Eq. (10) gives for the plate temperature 


Tag ; td Ve 
T = T, + 2[Pr/(2 — Pr)]' ‘(1 — =sin™'[Pr/2]' ‘) . (12) 
TT 
Thus for the infinite plate the recovery factor is 
¥ 71/2 ae 1/2 
r = 2(Pr/(2 — Pr)]**|1 — =sin™ [Pr/2] ‘), 
us 


which is to be compared to the corresponding recovery factor for a finite plate with a 
steady boundary layer 


r = Pr”, (13) 


For this latter result see for example Emmons and Bainerd’. Figure 1 compares these two 
results. The equilibrium temperature of an accelerated plate is independent of time and 
is very nearly the same as that of the plate in steady flow of the same velocity. The tem- 
perature distribution of Eq. (10) has been obtained by graphical integration and is shown 















in Fig. 2. 
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By comparing the thickness of the boundary layer of Eq. (2), and Fig. 2 with that of a 
point on a steadily moving plate distant x from the leading edge it is found that the 

*Emmons and Bainerd, Temperature effects in a laminar compressible fluid boundary layer along a flat 
plate, J. Appl. Mech., p. 4 (1941). 
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non-steady layer would equal the steady layer in times given approximately by 


5x ‘ 
( = 7 for velocity, 
62 
( = 7 for temperature. 


Thus to a reasonable approximation it can be said that by the time a point on a suddenly 
accelerated plate moves 5 times its distance from the leading edge, its boundary layers 
will have become steady state ones. 


PLASTIC WAVE PROPAGATION IN A BAR OF MATERIAL 
EXHIBITING A STRAIN RATE EFFECT" 


By L. E. MALVERN (Carnegie Institute of Technology) 


1. Introduction. The propagation of a transient wave of plastic deformation due to 
longitudinal impact on a bar has been treated by Donnell,* and White and Griffis,* by 
a non-linear superposition method. The partial differential equations governing the 
wave propagation were derived independently by Taylor* and von Karman’ under the 
assumption of a relation between stress and strain independent of strain rate. Constant 
velocity tension impact tests at the California Institute of Technology’ gave fair 
agreement with the theory. Some systematic discrepancies were, however, observed. 
In the tension impact tests the maximum residual strain was smaller than predicted 
by the theory, and the observed force-time variation at the fixed end during impact 
showed that the stress there was greater than the theory predicted. It has been sug- 
gested” that these discrepancies were due to the use in the theory of an invariant relation 
between stress and strain independent of strain rate. At the high strain rates involved 
in deformation under impact a considerable deviation from the static stress-strain 
relation may be expected. The present work extends the theory to apply to materials in 
which the stress is a function of the instantaneous plastic strain and strain rate. 


‘Received June 5, 1950. The results presented here were obtained in the course of research conducted 
at Brown University under Contract N7onr-358 sponsored jointly by the Office of Naval Research and the 
Bureau of Ships. This paper is part of a thesis submitted in partial fulfillment of the requirements for 
the degree of Doctor of Philosophy at Brown University, October, 1949. 

2L. H. Donnell, Longitudinal wave transmission and impact, A.S.M.E., Trans. 52 (1), APM 153-167 
(1930). 

8M. P. White and L. Griffis, The permanent strain in a uniform bar due to longitudinal impact, J. Appl. 
Mech., A.S.M.E., Trans. 69, A-337-A-343 (1947). 

‘G. I. Taylor, Propagation of earth waves from an explosion, British Official Report R.C. 70 (1940). 

‘Th. v. Karman, On the propagation of plastic deformation in solids, N.D.R.C. Report No. A-29 
(O.S.R.D. No. 365) (1942). 

6P. E. Duwez, D. 8. Wood, D. S. Clark, and J. V. Charyk, The effect of stopped impact and reflection 
on the propagation of plastic strain in tension, N.D.R.C. Report No. A-108, (O.S.R.D. No. 988) (1942). 

7P. E. Duwez and D. 8. Clark, An experimental study of the propagation of plastic deformation under 
conditions of longitudinal impact, A.S.T.M., Proc. 47, 502-532 (1947). 
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2. Proposed flow law for tension and compression. It is assumed that in longi- 
tudinal impact on a cylindrical or prismatic bar a relation of the form 

c= o(e’’, ée’’) (1) 

exists between the values of the nominal tensile stress o (longitudinal force per unit of 

initial cross-sectional area), plastic strain ¢’’ (permanent change in length per unit 

initial length), and the plastic strain rate ¢’’’. Since ¢ is in general an increasing function 

of «’”’ this determines e’” as a function of o and e”’. This relation may be expressed as 


E,e'”’ = g(a, €), 


where the factor E, is Youngs modulus and « is the total strain. Elastic deformation is 
assumed to be independent of strain rate. Thus, if e’ denotes the elastic strain, 


E,¢’ = oa’. (2) 
The relation between total strain, strain rate, and stress is then 
E,e = o° + g(a, €). (3) 
The static stress-strain relation ¢ = f(e) is interpreted as a succession of equilibrium 
states so that plastic flow occurs only when the plasticity condition 


a > fle) (4) 


is satisfied. Otherwise the elastic law (2) applies instead of the plastic flow law (3). 
The elastic law also applies until the initial yield strain ¢, is reached on the first loading. 
The plasticity condition (4) as stated applies to tensile impact (¢ and e positive) but 
the same form of law may be used in compressive impact if compressive stress and strain 
are reckoned positive. 

There is some evidence that the right-hand member of (1) should have the form of 
J(€) plus a term depending logarithmically on the plastic strain rate.’ If this is the case 
g(c, €) will depend exponentially on o — f(e), the excess of the instantaneous stress 
over the static stress at the same strain. Sokolovsky’ has treated wave propagation 
in a material without work-hardening using a law of the form (3) in which g(a — o,) 
was a function only of the excess of the stress over the initial yield stress a, . 

3. Equations governing the wave propagation. The propagation of the wave of 
plastic deformation is governed by the following system of three partial differential 
equations, in which zx denotes the initial distance of a cross section from the impact 
end, v the particle velocity (assumed constant over each cross section, and p the initial 


density. 
00 Ov 
ee ou *® 
de Ov 
a” (5) 
Oe eleg 
Eo at =" at = g(c, €). 


8See, for example, H. Deutler, Experimentelle Untersuchung ueber die Abhaengigkeit der Zugspannun- 
gen von der Verformungsgeschwindigkeit, Phys. Z. 33, 247-259 (1932). 
°V. V. Sokolovsky, The propagation of elastic-viscous-plastic waves in bars, Prikl. Mat. i Mek. 12, 


261-280 (1948). 
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The first equation is the longitudinal equation of motion. The second equation is a 
consequence of the fact that v = du/dt and e = du/dx, where u(z, t) is the displacement 
at time ¢ of the cross section which was initially at distance x from the impact end of 
the bar, and the third equation is the law (3). 

The system (5) is a hyperbolic system of quasi-linear partial differential equations 
which may be integrated numerically by the method of characteristics under appropriate 
boundary conditions. The characteristics in the z,t-plane are the three families of straight 
lines defined by the characteristic differential equations 


dx = 0, dz — c, di = 0, dx + ce dt = 0, (6) 
where 
Co = (E,/p)*”” 


is the constant speed of propagation of longitudinal elastic waves in the bar. 
The following three equations hold respectively along the three characteristics 
defined by (6). 


E,de—do= gc, dt 
do — pct, dv = — g(a, € dt (7) 


do + pt, dv = — gle, €) dt 


The equations (6) may be integrated immediately to give the fixed straight character- 
istics of the plastic region of the z,t-plane, but the equations (7) will in general require 








UNDISTURBED REGION 


X 





Fig. 1 Characteristics in the z,t-plane. 


step-by-step numerical integration. For this purpose the differentials of (7) are replaced 
by finite differences, and the value of g(c, «) in each equation is taken as an average of 
the values along the appropriate segment of the characteristic. 

4. Boundary conditions for continued impact on a semi-infinite bar. Consider a 
continued tensile impact on a semi-infinite bar which is initially at rest. The extension 
of the theory to finite bars and finite durations of impact is not difficult, although the 
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numerical integration becomes more involved. The z-axis is chosen so that at the instant 
of impact, the impact end is at the origin and the bar lies along the positive x ax.s. Suppose 
that at time ¢ = 0, the end of the rod is instantaneously set in motion with velocity 
v) such that | v | > coe, after which the end velocity is maintained constant at the value 
vo . For a tensile impact the velocity is negative. 

Immediately upon impact a shock wave of elastic deformation begins to travel 
along the bar at the speed c, . This leading wave front is represented in Fig. 1 by the 
line x = ct. The shock wave conditions which hold across an elastic shock wave traveling 
in the positive direction are 


Ac = —pc,Ar, 
Av = —c,Ae, (8) 


Ao = pcoAe = E, Ae, 


where Ac, Ae, and Av are the jumps in stress, strain, and velocity, respectively, as the 
shock wave passes. The first condition results from equating impulse to change of 
momentum for the traversing of an element of the bar by the shock wave. The second 
condition is a consequence of the continuity of displacement across the shock and the 
third condition follows from the first two. 

Since o = e€ = v = O in the undisturbed region ahead of the shock wave, the jump 


conditions (8) yield 


c= plo€ = — ply on r= Col (9) 
just after the shock wave passes. The line x = ef is also a characteristic of (5) along 
which the second equation of (7) holds. This equation may be integrated after elimi- 


nating ¢ and v by use of (9) to yield 


P* dr ! 
| ——*_,, = <1, (10) 
Ja, G\T, T/ pl « 
where o) = —pCo’o is the stress at x = 0, t = 0. Equations (8) and (10) thus determine 


o, ¢, and v along x = et . With this information and the boundary condition v = v, 


on x = 0, the numerical integration of the equations (7) may be performed to determine 
o, e and v throughout the plastic region. The values so obtained should be checked 


with the plasticity condition (4) at each point to make sure that the point is in the 


plastic region. 
5. Impact of finite duration: unloading. If the impact is of duration ¢, , after which 


the boundary condition is ¢ = 0 on x = O, the solution may be constructed as in the 
preceding section up to the characteristic of the family dx — c, dt = O passing through 
the point (0, ¢,), i.e. up to the line MN of Fig. 2. In the unloading region Eq. (2) replaces 
Eq. (3) so that g(c, «) = 0 for unloading. Equations (6, 7) then yield 

o — Eve = const. on x = const., 

o — pC) = const. on x — Cot = const., (11) 

o + pcv = const. on x +c t = const. 
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The various constants are determined by matching solutions at the elastic-plastic un- 
loading boundary in the z,t-plane. 


+ 











X 


0 


Fig. 2 Unloading shock wave MN after impact of duration t on a semi-infinite bar. 


A sudden reduction of the impact end stress to zero causes the initial unloading wave 
to be a shock wave. Since the shock wave is elastic the jump conditions (8) apply and 
the wave travels at the speed ¢y , i.e. along MN, Fig. 2. The shock intensity decreases 
as the wave progresses and may decrease to zero. If the shock wave travels as far as 
Q, Fig. 2, the solution is determined in the triangle MQR by Eqs. (11), the boundary 
condition ¢ = 0 on x = O, and the jump conditions (8) together with the known values 
of o, ¢, and v along MN before the shock wave passes. 

If the shock wave is absorbed, say at Q, the elastic-plastic boundary becomes a 
boundary of continuous transition from the plastic to the elastic state, and does not in 
general continue along the line /N. Along a continuous unloading boundary the static 
relation « = f(e) holds. This condition and the characteristic conditions in the elastic 
and plastic regions suffice to determine the position of the unloading boundary, although 
a trial and error procedure is usually needed to find the boundary points. A continuous 
transition from plastic to elastic behavior may begin even before the impact ends. The 
unloading shock wave then travels through an elastic region until it overtakes the 
continuous unloading boundary. 

6. Solutions of the equations. Numerical integration of the system (5) has been 
carried out for an idealized form of the law (3) in which g(o, «) = k[o — f(e)] where k 
is a multiplicative constant, with the function f(¢€) chosen in a simple form which ap- 
proximates the static stress-strain curve of a hardened aluminum. Even with the idealized 
law used the solutions indicate that a law of the type (3) can account for the discrepancies 
observed in the stress-time variation at the fixed end of impact specimens. This type 
of law does not, however, account for the discrepancy observed in the maximum residual 


strain. 
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It was possible to obtain an explicit solution for the equations in the case of a non 
work-hardening material of a type previously treated numerically by Sokolovsky.’ With 


g(c, €) = k(o — o,), where o, is the yield stress, and with the introduction of the di- 
mensionless variables 
T = ki X = rk/e 
S = a/o, E = ¢/e, V = v/ec¢, ; 
the system (5) becomes 
OS _ aV 
aX aT’ 
OE _ aV 
oT = aX’ 
dE as ; 
am = am S— I. 
oT ~ ar * 
Elimination of EZ and V yields 
oS oS a8 
a —  - =. 


ax’ aT” aT 
The boundary conditions for S in the plastic region of the X,T-plane (above X = T) are 


aS ‘ 
+7 = ( on A=Q9@ 
S = 1+ C exp (-—T/2) on X = T, 
where 
C=|/V_| — 1 = | 19/ece, | — 1. 
Introduction of the canonical variables, a = 3(7.+ X) and 8 = 3(T — X) and a new 
function P defined by 
S = 1+ CP exp[—}(a + 8)] 


yields the telegraph equation for P, 


ni png 
da 0B 4 
with boundary conditions 
oP 
= =0 on a=86 
on 


Pat on B=Q0O 


By a Laplace transform technique, or by the Riemann method of integration it may be 


determined that 
P(a, 8) = I,[(aB)"””] 
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is the solution required, where J, is the Bessel function of order zero with pure imaginary 
argument. Thus, returning to the variables, X and T, the dimensionless stress solution is 


S(X, T) = 1+ CI,[4(T? — X*)'”] exp (—T/2). 


The strain and velocity solutions are then obtained from the characteristic conditions 
by numerical quadratures. 


BOOK REVIEWS 


Nonlinear vibrations in mechanical and electrical systems. By J. J. Stoker, Interscience 
Publishers, Inc., New York, 1950. xv + 273. $5.00. 


Chapter I contains a concise review of linear vibrations in one degree of freedom. This is followed, in 
Chapter II, by a discussion of free vibrations of undamped systems with nonlinear restoring forces. One 
finds here, in addition to the material given in the books by Timoshenko (Vibrations Problems in Engin- 
eering. D. Van Nostrand Co., N. Y., 1937) and Den Hartog (Mechanical Vibrations. McGraw-Hill Book 
Co., N.Y., 1940) the notion of the energy curves in the phase plane. The theory of this chapter is illus- 
trated by interesting examples. 

Chapter III is divided into three sections; the first explains Liénhard’s graphical construction of 
trajectories in the phase plane, the second gives the theory of singular points of first order differential 
equations, and the third shows applications of that theory to nonlinear mechanical and electromechanical 
systems. A good deal of attention is devoted to a novel treatment of the elastic stability of columns from 
a dynamical viewpoint, making use of the notion of singular points. 

Chapter IV deals essentially with Duffing’s equation and its integration by both an iteration and a 
perturbation method. In this chapter reliance is placed on physical intuition, the more delicate questions 
of the existence of the perturbation series, and the stability of solutions being left for later treatment. A 
special case of subharmonic resonance in Duffing’s equation with and without damping is discussed, and 
Rauscher’s iteration method is briefly described. Following Rauscher’s method, combination tones are 
introduced, i.e., response frequencies which are linear combinations of the frequencies w; and w. of a forcing 
function H, cos wt + H: cos wt. The chapter ends with some introductory remarks regarding the stability 
of harmonic solutions and with a helpful table summarizing the essential differences in the response of 
Duffing’s equation and of the related linear equation. 

Chapter V treats, in two sections, Van der Pol’s equation without and with forcing function. In the 
first part the autonomous equation is discussed and a very appealing treatment of relaxation oscillations, 
and also a higher approximation are included. The second part deals with Van der Pol’s method applied 
to the non-autonomous system and with the stability investigations by Andronow and Witt. It gives 
altogether a very complete treatment for the two cases of a forcing frequency near the natural frequency 
of the system, and not near that frequency. 

Chapter VI, the last chapter of the book, treats linear equations with periodic coefficients because of 
their applicability to investigations of the stability of solutions to nonlinear problems. The Floquet theory 
is well reproduced as is an interesting and simple discussion of the stability of the Mathieu equation 
w’’ + (6 + «cos z)w = 0 in the neighborhood of « = 0. Finally, the stability of harmonic solutions of 
Duffing’s equation without damping is investigated by the theory of this chapter. 

This final chapter is followed by six appendices giving, in that order, (1) the mathematical justifica- 
tion of the perturbation method, (2) the existence of combination oscillations, (3) the existence of limit 
cycles in self-excited systems, (4) a proof of a previously introduced and attractive argument regarding 
the form of the limit cycle of dv/dt = &[F(v) — £]/v for a special form of F(v), and as e + ©@ (called the 
limit of limit cycles), (5) Poincaré’s criterion for orbital stability, and (6) a very palatable proof of unique- 
ness of a limit cycle in the free oscillations of a self-sustained system. The last chapter and the appendices 
are distinctly not elementary in character. 
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Last, but far from least, the appendices are followed by a very helpful bibliography which makes 
reference not only to the source material of the book, but also to many topics only covered sketchily or 
merely mentioned without description. 

The author has, on the one hand, purposely omitted certain important general methods and results 
achieved in the field of nonlinear vibrations (see the ‘‘Introduction’’) because they can be found else- 
where. On the other hand, he has included a number of refinements which, although interesting, do not 
add a great deal to an understanding of nonlinear vibration phenomena. The deletion of general theories 
from a text for the reason that they can be found elsewhere is not conventional in textbooks. The reviewer 
believes that this book on nonlinear vibrations would have gained from a presentation of the methods of 
Kryloff and Bogoliuboff (Princeton University Press, 1943), the classical work on stability by Liapounoff 
(Annales de Toulouse, Paris, 1907), and at least a statement of certain general and useful theorems proved 
by Levinson (Journal of Math. and Phys., June, December 1943) and by Levinson and Smith (Duke 
Math. Journal, 1942). 

For reasons unclear to the reviewer, and quite unessential as the author observes, he has rendered 
Duffing’s equation so that the amplitude of the forcing function is small with the coefficient of the non- 
linear term. This appears to reduce the generality of the results. 


Among the sections in the book more difficult to follow one might mention the one on subharmonic 


resonance. In searching for a solution with frequency w/n, n > 1 being an intezer, and w being the forcing 
frequency. the author has chosen without further explanation n = 3. The reader is left wondering 
whether there is something magic about n = 3, and whether the treatment is the same for n even and n 
odd. This section suffers also from a lack of correspondence between figures and attending text, and from 
the use of the term bifurcation without any explanation of its meaning. 

As may be expected in a first edition, there are a few errors and many misprints in the book. It is 
hoped that these will be corrected in later editions. 

These remarks are not exhaustive, nor are they intended to detract from the merits of Stoker’s 
“Nonlinear Vibrations’’. The reader will find this book an excellent addition to the growing literature on 
nonlinear phenomena. The material selected for presentation is most interesting, and some of it seems to 
have been published here for the first time. The style is very readable, and the absence of a great mass of 
unessential mathematical detail is a welcome feature. 


R. M. RosenBERG 


Matrix analysis of electric networks. By P. Le Corbeiller. Harvard Monographs in Applied 
Science, Number 1. Harvard University Press, Cambridge, and John Wiley & Sons, 
Inc., New York, 1950. xi + 112 pp. $3.00. 

This is the first of a new series of monographs published by Harvard University under the care of an 
Editorial Committee (chairman, F. V. Hunt). In a Foreword the Committee states the purpose of the 
series: it is to be “a medium for publishing the results of University research to a wider audience than 
would be reached by individual professional journals”, and will be concerned primarily with the applied 
physical sciences, with emphasis on “‘intellectua] borrowing among the academic disciplines.’ 

A book must be written primarily for a certain group of readers, and this book has been written, 
naturally enough, for electrical engineers. But for review it has fallen into the hands of an applied mathe- 
matician who certainly does not know all he ought to know about networks in order to appreciate the 
finesses, and who perhaps is a little too much inclined to try to dig down to mathematical essentials. 

The book is pleasantly compact, nicely printed with many diagrams, written in an easy graceful 
style, and supplied with an index. It is unified in purpose, namely, “. . . to offer to engineering circles a 
simple approach to Gabriel Kron’s method of analysis of stationary electric networks”’. 

Chapter I disposes of preliminary concepts: definitions of terms, elements of matrix algebra against an 
electrical background, and a page and a half on “‘The Two Fundamental Topologic Equations’’. Chapter 
II, on the mesh method of Kron, contains the kernel] of the book. Chapters III and IV, on the node-pair 
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method and the mixed method respectively, take the reader into more complicated technical situations 
in order to show the versatility of Kron’s method. 

Since his purpose is to explain a method which has, in Kron’s own exposition, been found difficult to 
understand, the author rightly avoids complexities and deals with d.c. and a.c. currents in steady states. 
This means that the solution of a network (i.e. the determination of the branch currents in terms of the 
branch emfs) is nothing more than the solution of a set of linear algebraic equations. 

That sounds very simple, and in a sense it is: anyone, given a network and a pair of basic laws, can 
get all the equations he needs, and anyone with sufficient time and patience can solve a set of linear alge- 
braic equations. Kron’s method systematizes this crude procedure, and probably the best way to see what 
the method consists of is to write down the fundamental relations: (1) i = Ci’, e’ = Cie; (2) e = Z'7’, 
Z’ = C,ZC; (3) t = C(C.ZC) Ce. Here 7 is a column matrix of branch currents (B of them), e a column 
matrix of branch emfs (B of them), 7’ a column matrix of mesh currents (M of them), e’ a column matrix 
of mesh emfs (7 of them), Z a B X B impedance matrix, and finally C (the central theme) is Kron’s 
rectangular transformation matrix and C;, its transpose. 

Suppressing for a moment our curiosity as to the source of these relations, let us look at the last (3). 
It is precisely what we want, a formula expressing the branch currents in terms of the branch emfs, and 
so solving our problem. All we need to know are (7) the impedance matrix Z and (77) the transformation 
matrix C. The rest is routine matrix computation. 

We do not have to bother about Z, for it is supposed given to us. It is C we must look into. Its mean- 
ing is demonstrated with admirable clarity on pp. 27-33 of the book in terms of a simple network with four 
nodes and six branches. Intuition and reason walk hand in hand in this simple illustration and the reader 
turns over from p. 33 to p. 34 with his curiosity put on edge by the promise that the fundamental relations 
will now be proved to hold for all cases. 

The problem is clear. In our equation (3) above the meanings of the symbols 7, e, Z are clear, and we 
are, I believe, clear as to what we mean by a network. Then, if vou draw a network for me and give me 
the square matrix Z, can I, or can I not, write down a matrix C such that (3) is the true connection be- 
tween branch currents and emfs? I accept the cha!lenge—draw me a network with a thousand nodes and 
ten thousand branches, and I shall find you a matrix C! 

Why am I so confident? When I look into the problem, I realise that it has really nothing to do with 
Z: the equations (1) above are all that matter, and the second is simply the first re-dished, so that atten- 
tion focuses on 7 = C7’. Once that is granted, all the rest follows simply in a few lines. And this equation is 


simply a mathematical statement to the effeet that in any given network it is possible to pick out a set of 


i 
meshes such that if arbitrary currents are assigned to them, then the branch currents obtained from these 
mesh currents by the ordinary method of superposition satisfy Kirchhoff’s first law at the junctions or 
nodes, and that in this way we can get the most general system of branch currents satisfying the law. The 
question is one of pure topology, covered by a transparent electric veil. 

My belief in the validity of i = C7?’ is primarily an act of faith. If it were not generally true, Maxwell’s 
mesh method would be nonsense, and Maxwell was no fool. And there are many electrical engineers, 


including Kron and Le Corbeiller, for whose intuition I have the highest respect; if there was something 


wrong here, they would have found it out. 
And vet ore cannot base theories on acts of faith; one asks for reasons and proofs. Having been too 
lazy to drill through the iron mask of topology to find the proof of i = C7’ which I know lies there, and 


having been even too lazy to digest properly the proof given by Ingram and Cramlet in the Journal of 
Mathematics and Physies, 23, 134-155 (1944), I turned over from p. 33 to p. 34 expecting to find at least a 
plausible proof of the fundamental formula on which Kron’s method rests. 
ut I did not find it. Instead I found the problem obscured with the primitive networks and inter- 

mediate networks of Kron, and I emerged with my faith unshaken but with my reason as virginal as when 
I started 

The truth of the matter seems to be that at the basis of network theory there lies an element of 
topology, perhaps a small matter to a professional topologist, but formidable to an electrical engineer or 
to an applied mathematician who is repelled by the formal exactitude of topology as it is currently treated. 
Professor Le Corbeiller disclaims any intention of proving topological formulae, and one can sympathise 
with that attitude. But it is an impossible attitude to adopt if one is to use the word “proof” at all in 
connection with what is an essentially topological subject 

However, ‘“‘proof”’ is a subjective word and its mathematical meaning is chilling to the rich intuitive 
understanding which is the essence of engineering. Rigid reasoning is not necessarily the same thing as 
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understanding, and it may well be that Professor Le Corbeiller has chosen the right path to his goal. As 
far as one reader at least is concerned, he has achieved it; I understand now what Kron’s method is, and 
I owe that understanding to Professor Le Corbeiller’s excellent exposition of it, even if I have to look 


elsewhere for the “‘proofs’”’. 
J. L. SYNGE 


Philosophy of mathematics and natural science. By Hermann Weyl. Princeton University 


Press, 1949. 371 pp. $5.00. 


In 1926 Mr. H. Wey! contributed to R. Oldenbourg’s Handbuch der Philosophie a substantial article 
in two parts dealing, the first one with the philosophy of mathematics, the second one with that of physics. 
The present volume is a translation of this article, brought up to date at certain points, and enriched 
with six essays which represent the author’s present thought in a field now extending into chemistry and 


biology. 

This is a book no one should miss who likes science and philosophy. I do not say ‘philosophy of 
science”’, a fashionable topic at the present time. This is a book by a master of mathematical physics who 
also knows philosophy, classical philosophy, philosophy proper. To read and to reread it is a profit and a 


delight. 

There is a marvelous bibliography, or rather, which is even better, a dozen lists or so, at irregular 
intervals, of books the author has liked and which he has quoted, or might have quoted, in the last section 
or chapter. This is an invitation to browse in Mr. Weyl’s own study of which every thoughtful reader 
will be duly appreciative 


And the index is excellent. 
P. LECoRBEILLER 


Analysis and design of experiments. By H. B. Mann. Dover Publications Inc., New York, 
1949. 195 pp. $2.95. 


An interesting aspect of the development of statistics during the past quarter-century is the growth 
of a field usually known as the ‘design of experiments.” In the larger centers of statistical teaching, a 
course in this topic is now a standard part of the curriculum. A word of explanation as to how this arose 
is, I think, necessary to an understanding of the objectives of Dr. Mann’s book. 

To quote a description given by Yates, the essence of an experiment is the “imposition of deliberate 
change with the intention of studying its effects.’’ In most types of experimentation, however, we cannot 
measure exactly the effects of an imoosed change, because these effects are influenced by extraneous 
variations, to which the term ‘‘experimental errors”’ is often given. The presence of extraneous influences 
of this type greatly complicates the problem of drawing sound conclusions from a completed experiment, 
for the results are found to vary from one part of the experiment to another, so that it may not be clear 
just what can be inferred with confidence. In this situation experimenters have often consulted statis- 
ticians, because of the statistician’s knowledge of the theory of errors and of his preoccupation with the 
drawing of inferences from variable data. That is, the statistician was being consulted about the analysis 


of the results, not about the conduct of the experiment. 

Needless to say, the way in which an experiment is conducted determines the type of conclusion, if 
any, that can be drawn from its results. Consequently the statistician was led naturally to investigate 
the logical principles that must govern the conduct of an experiment if it is to permit reasonably clear 
inferences. This study was initiated about 25 years ago by R. A. Fisher and its fruits form the basis of 
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his book ‘“‘The Design of Experiments’’. Two principal developments may be mentioned. The first is an 
elegant method of computation known as the analysis of variance, by which the calculations needed for 
performing tests of significance and constructing confidence limits are provided. The analysis of variance 
is a particular application of the method of least squares, but it can be learned and used competently by 
experimenters who would disclaim any knowledge of probability or mathematics. The second outcome 
is a series of ingenious methods for carrying out experiments, discovered mainly by Fisher and Yates, 
which allow us to eliminate the influence of many extraneous sources of variation from the effects under 
investigation, with a resultant gain in the precision of the experiment. 

These two developments form the content of Dr. Mann’s book. As he points out, the analysis of 
variance is usually presented in the literature for digestion by experimenters rather than by mathe- 
maticians, so that the assumptions required and the mathematical justification of the theorems utilized 
are not fully set forth. The first six chapters give a lucid and adequate exposition of the theory of the 
analysis of variance as it applies in the common types of experimental design, starting with a knowledge 
only of calculus and the rudiments of probability. Later chapters present the extensions to nonorthogonal 
data, including that which arises in incomplete block designs, to factorial experiments, and to the analysis 
of covariance. 

The remainder of the book deals with certain problems, of interest to mathematicians, that have 
arisen in the construction of experimental designs. Complete sets of (m — 1) orthogonal latin squares of 
side m have proved useful. Dr. Mann gives a summary of the present state of knowledge on this subject. 
By use of the properties of Galois fields, he proves the well-known result that a complete orthogonal set 
can be constructed for any square whose side is the power of a prime. Beyond this, little is known except 
two results of restricted scope. The first is that if the side of the square is of the form IIp,;”", where the px 
are primes, we can construct a set of r orthogonal squares, where r = min. (p;"" — 1). The second result, 
on the negative side, is that we cannot construct a set of two orthogonal squares of side 6. 

Some account is given of the canstruction of balanced incomplete block designs. These are arrange- 
ments of ¢ different letters in blocks or groups of size k (k < t), such that every pair of letters occurs 
times within the same block. For experimentation, interest centers on arrangements in which ) is as 
small as possible. Dr. Mann shows how finite geometries have been applied in the construction of these 
arrangements. 

The book is designed for three purposes: for mathematicians with no background in statistics who 
wish to become acquainted with the theory underlying the design of experiments, for experimenters with 
the same objective, and for use as a text-book in a course on the subject. The exposition is admirably 
clear throughout. In a book of limited length, such as this, the topics discussed are necessarily restricted. 
For the benefit of mathematical readers, some examples of experiments in different fields of research 
might have been presented in detail, in order to bring out clearly why the analysis of variance is used in 
the summary of results, and why such arrangements as latin squares and balanced incomplete blocks help 
to solve the experimenter’s problems. As it is, such readers are likely to have a number of unanswered 
questions after a perusal of the book. For experimenters with the requisite mathematical background, the 
book should prove both useful and stimulating. 

W. G. CocHRran 


Practical analysis. By Fr. A. Willers. Translated by R. T. Beyer. Dover Publications, 
Inc., New York, 1948. 422 pp. $6.00. 


This German textbook covers very thoroughly the various types of numerical problems which the 
individual scientist or engineer is apt to encounter in the course of his work. Chapter 1 describes several 
types of tools such as slide rules and nomograms; the section on calculating machines was entirely re-writ- 
ten for the American edition by T. W. Simpson. Chapter 2 covers Interpolation; Chapter 3, Approximate 
Integration and Differentiation; Chapter 4, Algebraic Equations; Chapter 5, Analysis of Empirical Func- 
tions; Chapter 6, Approximate Integration of Differential Equations. 

A feature of this book is that almost as much space is given to graphical methods as to digital com- 
putation. This makes the treatment very pleasant and easy to follow. It is all the more regrettable that 
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the numerous and carefully executed diagrams of the German edition should have been, on the whole, 
very poorly reproduced. 

A number of mathematicians are at present servicing very high-powered electronic machines and the 
kind of problems in numerical analysis which they encounter is not considered in the present book. 
Treatises on numerical analysis which we hear are in preparation on several sides will probably concen- 
trate on those, to the benefit of the professional, and leave aside the individual worker. This one will still 
have to execute for himself, at more or less long intervals, some tabulation or integration work with the 
help of a modest desk-machine; he will find the book by Dr. Willers, with its enormous wealth of examples, 


a most helpful instructor and guide. 
P. LECORBEILLER 


Classical mechanics. By Herbert Goldstein. Addison-Wesley Press, Inc., Cambridge 42, 
Massachusetts, 1950. xii + 399 pp. $6.50. 


In this book the author succeeds in a very satisfactory way in showing how classical mechanics leads 
into the various branches of physics as well as in discussing the mathematical methods. 

The chapter headings are: 1. Survey of Elementary Principles. 2. Variational Principles and La- 
granges Equations. 3. The Two-Body Central Force Problem. 4. The Kinematics of Rigid Body Motion. 
5. The Rigid Body Equations of Motion. 6. Special Relativity in Classical Mechanics. 7. The Hamilton 
Equations of Motion. 8. Canonical Transformations. 9. Hamilton Jacobi Theory. 10. Small Oscillations. 
11. Introduction to the Lagrangian and Hamilton Formulations for Continuous Systems and Fields. 

The discussion of canonical transformations, integral invariants, Lagrange and Poisson Brackets, 
and the Hamilton Jacobi theory are treated in detail; the connection with geometrical optics and wave 
mechanics is concisely shown 

The discussion of central field problems includes some discussion of collision and scattering processes. 

The section on Kinematics of Rigid Body Motion makes use of matrix methods and treats in detail 
the use of the Cayley-Klein parameters for rigid body rotation; the connection with Pauli-Spin matrices 
is shown. 

The reviewer feels that this book will have a strong appeal to physicists, to physical chemists, and to 


some applied mathematicians. 
Roun TRUELL 


Wave guides. By H. R. L. Lamont. Methuen & Co., Ltd., London, and John Wiley & 
Sons, Inc., New York. vii + 118 pp. $1.50. 


This is the third edition of a Methuen monograph revised in 1950 and first published in 1942. A 
concise presentation of the essential parts of wave guide theory. It contains a discussion of wave guides 
as transmission lines, as resonators and as radiators. The book contains a brief outline of the diffraction 
problem in an aperture together with the use of W. R. Smythe’s results. There is a somewhat hasty men- 
tion of diaphragms in rectangular guides, change in guide cross section, junctions, slots, and corrugated 


guides. 
Roun TRUELI 


Electromagnetic waves. By F. W. G. White. Fourth edition (revised 1950). Methuen & 
Co. Ltd., London, and John Wiley & Sons, Inc., New York. viii + 108 pp. $1.25. 


The first half of this Methuen monograph is assigned to a review of electromagnetic theory with a 
discussion of plane wave solutions of the field equations. The Lorentz dispersion theory is reviewed. The 
second half of the monograph is devoted to application of the dispersion theory to obtain the approxi- 
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mate optical properties of metals, the diffraction of x-rays, and the group velocity in a dispersive medium. 
The propagation in a dispersive medium with an applied magnetic field and the propagation of radio 
waves in the upper atmosphere are the topics of the last two chapters. Unfortunately the Bibliography on 


radio wave propagation contains nothing more recent than 1933 
Roun TRUELL 


Heaviside’s electric circuit theory. By H. J. Josephs. John Wiley & Sons, Inc., New York 
and Methuen & Co. Ltd., London. viii + 115 pp. $1.25. 


This is the second edition of a monograph first published in 1946. It is a very concise and well written 
exposition on the methods of Heaviside. It includes a rather complete discussion of the expansion theorem, 
operational formulae, and methods leading to integral equations with solutions by operational methods. 
Borel’s theorem, “‘shifting’’, and fractional derivatives are reviewed. Impulse functions are dealt with in a 
chapter entitled ““Heaviside’s Last Theorem’’. The application of these methods to transmission lines 
includes a discussion of the divergent expansions that arise. There is a final chapter on the application of 
contour integration methods to the solution of circuit problems. 

Roun TRUELL 


Radio aerials. By E. B. Moullin. The Clarendon Press, Oxford, 1949. xi + 514 pp. $8.00. 


This is the second volume in the series, “International Monographs on Radio,” edited by Sir Edward 
Appleton and Henry G. Booker. Times have changed! A few years ago not a single book could be found 
devoted entirely to antennas; now there are several excellent books, all coming from Great Britain, which 
treat this subject broadly, and some, such as the present volume, are concentrating on a single special 
class of antennas. There may be a difference of opinion as to whether antenna reflectors made up of planes 
and circular cylinders are of sufficient practical importance to deserve five hundred pages; but there is no 
question that this book contains much useful information, both theoretical and experimental. 

As far as theory is concerned, most problems considered in the book are two-dimensional and highly 
idealized. The author regards his theoretical solutions as a guide to experimental work and not as predic- 
tions of performance of practical reflectors. He chooses those problems which can be solved exactly and 
counts on the common sense of the engineer to derive the maximum benefit in planning experiments and 
anticipating the results of experiments. The experimental section takes one-third of the book and contains 
information on parabolic reflectors in addition to those considered in the theoretical section. 

The field of two-dimensional problems admitting exact solutions is covered exhaustively. The author 
considers infinitely long wires, systems of such wires, infinite gratings, infinite planes, wedges made with 
semi-infinite planes, wedges with cylindrical backs. The analysis of semi-infinite planes permits him to 
estimate the edge effects in reflectors of finite size. Of the three-dimensional problems only a few are dealt 
with: a dipole, a half-wave antenna, arrays of half-wave antennas, a small loop, a large loop carrying 
uniform current, and a circular ribbon-of dipoles. Radiation patterns constitute the principal subject of 
study. 

The author begins the book with a brief discussion of Maxwell’s equations, introduces retarded poten- 
tials, and obtains the wave equation. Then, he asks the reader to accept the solution of this equation and 
from it proceeds to develop the solutions for infinite current filaments by integration. Here the author has 
missed a unique opportunity to avoid the gap in the theoretical development without sacrificing the 
essentially elementary method of analysis. In the two-dimensional case, retarded potentials do not sim- 
plify solutions of problems; on the contrary, they complicate them. It is easy to obtain directly from 
Maxwell’s equations the expressions for the component of E parallel to the current flow—which the 
author does—and then relate these expressions to the given currents without bringing in the retarded 
potentials—which the author does not do. This approach would have enabled him to dispense with some 
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tiresome integrations, to avoid difficulties with divergent integrals, and to extricate himself from diffi- 
culties in connection with the dihedral reflectors whose angles are not submultiples of 360°. 

It is to be regretted that the author is prejudiced against what he terms the “optical approach’’ to 
problems of radiation. He admits that this approach is not wrong. Fairly enough he warns the reader 
that the book is written in the electrical language and that he will insist “on trying to relate fields to 
currents flowing in conductors.’’ In this he follows a well-established tradition in radio engineering and 
preserves a consistent point of view throughout the book. But the price of consistency seems high. Very 
useful aids to thinking and powerful methods for solving radiation problems are thus neglected. An 
opportunity to compare different approaches is lost. And the time is ripe for future radio engineers to 
learn the mysteries of languages other than the “electrical language.”’ 

This book was intended primarily for reference. Its subject is probably too limited to recommend it 
as a sole textbook; but students should find the book very instructive in that it shows by easy steps how to 


construct the solutions of various problems from simple basic solutions. 
S. A. SCHELKUNOFF 


Physics in Industry. The acceleration of particles to high energies. Based on a session ar- 
ranged by the Electronics Group at the Institute of Physics Convention in May, 1949. 
Institute of Physics, London, 1950. x + 58 pp. $1.60. 


This concise little book arose as a consequence of a convention at Buxton, England in May 1949. 
The book should be of real value to the non-specialist who is interested in the progress in the field of high 
energy accelerators, but who for some reason does not follow the literature. There are five chapters 
entitled: The Cyclotron (and synchro-cyclotron) (13 pages); Betatrons and Synchrotrons (16 pages); 
Electrostatic Generators (8 pages); Linear Accelerators (6 pages); Discussion of the Preceding Topics 
(14 pages). 

In general each section discusses briefly principles, design, and performance. Naturally, the dis- 
cussion of actual design and construction centers on the British work although the bibliographies seem 
to be rather complete. In fact, the references and bibliographies alone seem to the reviewer to be well 


worth the price of the book. 
Roun TRUELL 


Theory of oscillations. By A. A. Andronow and C. E. Chaikin. English language edition 
edited under the direction of Solomon Lefschetz. Princeton University Press, New 
Jersey, 1949. ix + 358 pp. $6.00. 


This volume is a translation and condensation of a much larger one in Russian. The reader should be 
thankful to Professor Lefschetz for his work of adaptation, pruning, and specific as well as general clarifi- 
cation. The work as it stands reads very easily. It deals mostly with autonomous (or isolated) oscillating 
systems of the flip-flop type, including the effect of the variation of one of the parameters of the system. 
Continuous systems are treated in the last chapter; non-linear systems subjected to applied forces, on 
which information is more easily available, receive passing mention. The reader is carried very pro- 
gressively through a series of examples, mechanical and electrical. These will take the conscientious 
reader a great deal of time to go through; however, this study should be very rewarding as it will allow 
him not only to understand fully the theoretical facts, but to build up his own philosophy of the subject. 

A group of ‘“‘advanced”’ fields: non-linear oscillations, feedback theory, Laplace and Fourier trans- 
forms, noise theory, and a few others, present quite a problem to the instructor and to the graduate 
student. On the one hand, it is impossible, without a good understanding of them, to analyze some of the 
most common and important phenomena occurring in practice. On the other hand, once the theory is 
well understood, its use is in many cases mostly qualitative or else requires much simpler calculations 
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than have been used while learning it. One would like to be able to impart to the students the essentials 
of these theories ‘‘in a nutshell’’, but personal experience cannot be by-passed. 

The solution would seem to be to encourage the better students to self-study of such subjects. The 
present book by Andronow, Chaikin (and Lefschetz) should be especially recommended for that purpose. 
The numerous examples treated have been collected in two separate tables at the end of the book, one for 
mechanical and one for electrical systems, thus making it possible for the student to restrict his introduc- 
tory work to the type of systems he is most familiar with. 

While the specialist will have his quarrels with the original authors about this or that minor point, 
the book as a whole should be commended as a most helpful contribution to a wider knowledge of non- 
linear systems, 

P. LeCorBEILLeR 























SUGGESTIONS CONCERNING THE PREPARATION OF 
MANUSCRIPTS FOR THE QUARTERLY OF 
APPLIED MATHEMATICS 


The editors will appreciate the authors’ cooperation in taking note of the following directions for the preparation of 
nanuscripts. These directions have been drawn up with a view toward eliminating unnecessary correspondence, 
avoiding the return of papers for changes, and reducing the charges made for “‘author’s corrections.” 





Manuscripts: Papers should be submitted in original typewriting on one side only of white paper sheets and be 
double or triple spaced with wide margins. Marginal instructions to the printer should be written in pencil to 
stinguish them cle: irly from the body of the text. 
The papers should be submitted in final form. Only typographical errors may be corrected in proofs; composi- 
tion charges for all major deviations from the manuscript will be passed on to the author. 


Titles: The title should be brief but express adequately the subject of the paper. The name and initials of the 
author should be written as he prefers; all titles and degrees or honors will be omitted. The name of the organization 
with which the author is associate d should be given in a separate line to follow his name. 

Mathematical Work: As far as possible, formulas should be typewritten; Greek letters and other symbols not 
avuilable on the typewriter should be ec arefully inserted in ink. Manuscripts containing pencilled material other 
than marginal instructions to the printer will not be accepted. 

The difference between capital and lower-case letters should be clearly shown; care should be taken to avoid 
confusion between zero (0) and the letter O, between the numeral one (1) ,the letter / and the prime (’), between 
alpha and a, kappa and &, mu amd uw, nu and 2, eta and n. 

The level of subscripts, exponents, subscripts to subscripts and exponents in exponents should be clearly 
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